#include <stdio.h>
#include <algorithm>
#include <vector>
using namespace std;

const long long mod = 1000000007;
const long long inv2 = (mod + 1) / 2;
vector<int> G[2020];

vector<long long> mul(vector<long long> &a, vector<long long> &b)
{
	vector<long long> r(a.size()+b.size()-1);
	for (int i=0;i<a.size();i++) for (int j=0;j<b.size();j++){
		r[i+j] = (r[i+j] + a[i] * b[j]) % mod;
	}
	return move(r);
}

vector<long long> add(vector<long long> &a, vector<long long> &b, int sig = 1)
{
	vector<long long> r(max(a.size(),b.size()));
	for (int i=0;i<a.size();i++) r[i] = a[i];
	for (int i=0;i<b.size();i++) r[i] = (r[i] + b[i] * sig % mod + mod) % mod;
	return move(r);
}

vector<vector<long long> > go(int x)
{
	vector<vector<long long> > r(3,vector<long long>(2,0));
	r[0][1] = 1;

	for (auto &y : G[x]){
		auto q = go(y);
		auto s = q[0];
		s = add(s,q[1]);
		s = add(s,q[2]);

		vector<vector<long long> > nr(3);
		for (int k=0;k<3;k++) nr[k] = mul(r[k],s);
		long long coeff[2][2] = {2,1,1,inv2};
		for (int i=0;i<2;i++) for (int j=0;j<2;j++){
			auto t = mul(r[i],q[j]);
			t.erase(t.begin());
			nr[i+1] = add(nr[i+1],t,-coeff[i][j]);
		}
		r = move(nr);
	}

	return move(r);
}

int main()
{
	int N; scanf ("%d",&N);
	for (int i=2;i<=N;i++){
		int p; scanf ("%d",&p);
		G[p].push_back(i);
	}

	auto r = go(1);

	long long ans = 0, f = 1;
	for (int i=1;i<=N;i++){
		f = f * i % mod;
		ans = (ans + (r[0][i] + r[1][i] + r[2][i]) * f) % mod;
	}
	printf ("%lld\n",ans);

	return 0;
}