#include <bits/stdc++.h>
using namespace std;
const int MAXN = 2010;
const long long MOD = 1000000007;
inline int read() {
  int x = 0, f = 1;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) x = (x * 10) + (ch ^ 48);
  return x * f;
}
inline long long qpow(long long a, long long b) {
  long long res = 1;
  while (b) {
    if (b & 1LL) res = res * a % MOD;
    b >>= 1, a = a * a % MOD;
  }
  return res;
}
int e, st[MAXN], to[MAXN];
int nxt[MAXN];
inline void Add(int u, int v) {
  to[++e] = v, nxt[e] = st[u];
  st[u] = e;
}
int n, fa[MAXN], size[MAXN];
long long dp[MAXN][MAXN][3], inv2;
long long tmp[MAXN][3];
long long fac[MAXN], ans;
inline void update(long long &cur, long long val) {
  cur = cur + val;
  if (cur >= MOD) cur -= MOD;
}
void dfs(int u) {
  int i, j, k;
  dp[u][1][0] = 1, size[u] = 1;
  for (i = st[u]; i; i = nxt[i]) {
    int v = to[i];
    dfs(v), size[u] += size[v];
  }
  if (size[u] == 1) return;
  for (i = st[u]; i; i = nxt[i]) {
    for (k = size[u]; k >= 1; k--) {
      tmp[k][0] = dp[u][k][0];
      tmp[k][1] = dp[u][k][1];
      tmp[k][2] = dp[u][k][2];
      dp[u][k][0] = dp[u][k][1] = dp[u][k][2] = 0;
    }
    int v = to[i];
    for (k = size[u]; k >= 1; k--) {
      for (j = min(size[v], size[u] - k + 1); j >= 1; j--) {
        if (j + k - 1 != size[u]) {
          update(dp[u][j + k][0], dp[v][j][0] * tmp[k][0] % MOD);
          update(dp[u][j + k][1], dp[v][j][0] * tmp[k][1] % MOD);
          update(dp[u][j + k][2], dp[v][j][0] * tmp[k][2] % MOD);
          update(dp[u][j + k][0], dp[v][j][1] * tmp[k][0] % MOD);
          update(dp[u][j + k][1], dp[v][j][1] * tmp[k][1] % MOD);
          update(dp[u][j + k][2], dp[v][j][1] * tmp[k][2] % MOD);
          update(dp[u][j + k][0], dp[v][j][2] * tmp[k][0] % MOD);
          update(dp[u][j + k][1], dp[v][j][2] * tmp[k][1] % MOD);
          update(dp[u][j + k][2], dp[v][j][2] * tmp[k][2] % MOD);
        }
        update(dp[u][j + k - 1][2], dp[v][j][0] * tmp[k][1] % MOD);
        update(dp[u][j + k - 1][2], dp[v][j][1] * tmp[k][1] % MOD * inv2 % MOD);
        update(dp[u][j + k - 1][1], dp[v][j][0] * tmp[k][0] % MOD * 2 % MOD);
        update(dp[u][j + k - 1][1], dp[v][j][1] * tmp[k][0] % MOD);
      }
    }
  }
}
int main() {
  int i;
  n = read();
  for (i = 2; i <= n; i++) Add(read(), i);
  inv2 = qpow(2, MOD - 2);
  dfs(1);
  fac[0] = 1;
  for (i = 1; i <= n; i++) fac[i] = fac[i - 1] * i % MOD;
  for (i = 0; i < n; i++) {
    int f = i & 1 ? MOD - 1 : 1;
    update(ans,
           f * fac[n - i] % MOD *
               ((dp[1][n - i][0] + dp[1][n - i][1] + dp[1][n - i][2]) % MOD) %
               MOD);
  }
  printf("%lld\n", ans);
  return 0;
}
