#include <bits/stdc++.h>
using namespace std;

const int N = 2005;
const int mod = 1e9 + 7;

int n;
vector <int> G[N];
int dp[N][N][3];
int f[N][3];
int nchild[N];

void add(int &x, int y) {
	x += y; while(x >= y) x -= mod; while(x < 0) x += mod;
}

void dfs(int u) {
	dp[u][0][0] = 1;

	for (int &v : G[u]) {
		dfs(v);
		memset(f, 0, sizeof f);

		for (int ncut_u = 0; ncut_u <= nchild[u]; ++ncut_u) {
			for (int deg_u = 0; deg_u <= 2; ++deg_u) {
				for (int ncut_v = 0; ncut_v <= nchild[v]; ++ncut_v) {
					for (int deg_v = 0; deg_v <= 2; ++deg_v) {
						if (!dp[u][ncut_u][deg_u] || !dp[v][ncut_v][deg_v]) continue;

						int mul = 1LL * dp[u][ncut_u][deg_u] * dp[v][ncut_v][deg_v] % mod;

						// cut
						add(f[ncut_u + ncut_v + 1][deg_u], (deg_v == 1 ? 2LL : 1LL) * mul % mod);

						// concatenate
						if (deg_u <= 1 && deg_v <= 1) {
							add(f[ncut_u + ncut_v][deg_u + 1], mul);
						}
					}
				}
			}
		}

		nchild[u] += nchild[v];
		for (int ncut_u = 0; ncut_u <= nchild[u]; ++ncut_u) {
			for (int deg_u = 0; deg_u <= 2; ++deg_u) {
				dp[u][ncut_u][deg_u] = f[ncut_u][deg_u]; // update the values of dp[]
			}
		}
	}

	++nchild[u];
	for (int ncut_u = 0; ncut_u < nchild[u]; ++ncut_u) {
		for (int deg_u = 0; deg_u <= 2; ++deg_u) {
			if (deg_u == 2) {
				add(dp[u][ncut_u][deg_u], dp[u][ncut_u][deg_u]); // x 2
			}
			if (u == 1 && deg_u == 1) {
				add(dp[u][ncut_u][deg_u], dp[u][ncut_u][deg_u]); // x 2	
			}
		}
	}
}

int main() {
	ios_base::sync_with_stdio(false); cin.tie(0);
	cin >> n;
	for (int i = 2; i <= n; ++i) {
		int p; cin >> p; G[p].push_back(i);
	}

	dfs(1);

	int fact = 1;
	int res = 0;
	for (int ncut = 0; ncut <= n - 1; ++ncut) {
		int npair = n - 1 - ncut;
		int cur = 0;
		for (int deg = 0; deg <= 2; ++deg) add(cur, dp[1][ncut][deg]);

		fact = 1LL * fact * (ncut + 1) % mod; // (ncut + 1)!
		cur = 1LL * cur * fact % mod;

		if (npair % 2 == 0) add(res, cur); else add(res, -cur);
	}

	cout << res << endl;
}