#include <bits/stdc++.h>

using namespace std;

const int Mod = 1e9 + 7;
const int Nmax = 2005;
typedef long long ll;

int n, i, j, x;
vector<int> v[Nmax];
int dp[Nmax][Nmax][3], fact[Nmax], D[Nmax][3], w[Nmax];

void add_to(int &x, int y) { x+=y; if(x>=Mod) x-=Mod; }
int inm(int x, int y) { return (ll)x*y%Mod; }

void Merge(int node, int son)
{
    int i, j, dg, dgr;
    for(i=0; i <= w[node] + w[son]; ++i)
        for(j=0; j<3; ++j)
            D[i][j] = 0;

    for(i=0; i<=w[node]; ++i)
        for(j=0; j<=w[son]; ++j)
            for(dg=0; dg<3; ++dg)
                for(dgr=0; dgr<3; ++dgr)
                    if(dp[node][i][dg] && dp[son][j][dgr])
                    {
                        add_to( D[i+j][dg], inm(dp[node][i][dg], dp[son][j][dgr] * (dgr ? 2 : 1)) );

                        if(dg < 2 && dgr < 2)
                            add_to( D[i+j-1][dg+1], inm(dp[node][i][dg], dp[son][j][dgr]) );
                    }

    w[node] += w[son];
    for(i=0; i<=w[node]; ++i)
        for(j=0; j<3; ++j)
            dp[node][i][j] = D[i][j];
}

void solve(int node)
{
    w[node] = 1; dp[node][1][0] = 1;
    for(auto son : v[node])
    {
        solve(son);
        Merge(node, son);
    }
}

int main()
{
 //   freopen("input", "r", stdin);
  //  freopen("output", "w", stdout);
    cin.sync_with_stdio(false);

    cin >> n;
    for(i=2; i<=n; ++i)
    {
        cin >> x;
        v[x].push_back(i);
    }

    fact[0] = 1;
    for(i=1; i<=n; ++i) fact[i] = inm(fact[i-1], i);

    solve(1);

    int ans = 0;
    for(i=0; i<=n; ++i)
        for(j=0; j<3; ++j)
            if( (n-i) % 2 == 0 ) add_to(ans, inm( dp[1][i][j] * (j ? 2 : 1), fact[i] ) );
                else add_to(ans, Mod - inm( dp[1][i][j] * (j ? 2 : 1), fact[i] ) );

    cout << ans << '\n';
    return 0;
}
