#include <bits/stdc++.h>
using namespace std;
const int MAX_N = 2005;
vector<int> G[MAX_N];
int dp[MAX_N][MAX_N][3][2];
int st[MAX_N], fac[MAX_N], inv2;
void make() {
  fac[0] = fac[1] = 1;
  for (int i = 2; i < MAX_N; i++) {
    fac[i] = (long long)fac[i - 1] * i % 1000000007;
  }
}
inline int mod_pow(int a, long long b) {
  int res = 1;
  while (b) {
    if (b & 1) {
      res = (long long)res * a % 1000000007;
    }
    a = (long long)a * a % 1000000007;
    b >>= 1;
  }
  return res;
}
inline int add(int x, int y) { return (x + y) % 1000000007; }
inline int sub(int x, int y) { return (x + 1000000007 - y) % 1000000007; }
inline int mul(int x, int y) { return (long long)x * y % 1000000007; }
void dfs(int u) {
  st[u] = 1;
  dp[u][1][2][0] = 1;
  for (auto&(v) : (G[u])) {
    dfs(v);
    for (int i = st[u]; i >= 1; i--) {
      for (int j = st[v]; j >= 1; j--) {
        dp[u][i + j][2][1] = add(
            dp[u][i + j][2][1],
            mul(dp[u][i][2][0],
                (long long)(dp[v][j][0][0] + dp[v][j][1][0] + dp[v][j][2][0]) %
                    1000000007));
        dp[u][i + j][1][1] = add(
            dp[u][i + j][1][1],
            mul(dp[u][i][1][0],
                (long long)(dp[v][j][0][0] + dp[v][j][1][0] + dp[v][j][2][0]) %
                    1000000007));
        dp[u][i + j][0][1] = add(
            dp[u][i + j][0][1],
            mul(dp[u][i][0][0],
                (long long)(dp[v][j][0][0] + dp[v][j][1][0] + dp[v][j][2][0]) %
                    1000000007));
        dp[u][i + j - 1][1][1] =
            add(dp[u][i + j - 1][1][1],
                mul(mul(dp[u][i][2][0], dp[v][j][2][0]), 2));
        dp[u][i + j - 1][1][1] =
            add(dp[u][i + j - 1][1][1], mul(dp[u][i][2][0], dp[v][j][1][0]));
        dp[u][i + j - 1][0][1] =
            add(dp[u][i + j - 1][0][1], mul(dp[u][i][1][0], dp[v][j][2][0]));
        dp[u][i + j - 1][0][1] =
            add(dp[u][i + j - 1][0][1],
                mul(mul(dp[u][i][1][0], dp[v][j][1][0]), inv2));
      }
    }
    st[u] += st[v];
    for (int i = (int)(1); i < (int)(st[u] + 1); ++i) {
      for (int j = 0; j < (int)(3); ++j) {
        dp[u][i][j][0] = dp[u][i][j][1];
        dp[u][i][j][1] = 0;
      }
    }
  }
}
int main() {
  cin.tie(0);
  ios::sync_with_stdio(false);
  int n;
  cin >> n;
  for (int i = 0; i < (int)(n - 1); ++i) {
    int p;
    cin >> p;
    G[p - 1].push_back(i + 1);
  }
  make();
  inv2 = mod_pow(2, 1000000007 - 2);
  dfs(0);
  int ans = 0;
  for (int i = (int)(1); i < (int)(n + 1); ++i) {
    if ((n - i) % 2) {
      ans = sub(ans, mul((long long)(dp[0][i][0][0] + dp[0][i][1][0] +
                                     dp[0][i][2][0]) %
                             1000000007,
                         fac[i]));
    } else {
      ans = add(ans, mul((long long)(dp[0][i][0][0] + dp[0][i][1][0] +
                                     dp[0][i][2][0]) %
                             1000000007,
                         fac[i]));
    }
  }
  cout << ans << "\n";
  return 0;
}
