use std::io;

const MOD: i64 = 1_000_000_007;

fn mod_pow(mut base: i64, mut exp: i64, modu: i64) -> i64 {
    let mut result = 1;
    while exp > 0 {
        if exp % 2 == 1 {
            result = result * base % modu;
        }
        base = base * base % modu;
        exp /= 2;
    }
    result
}

fn main() {
    let mut input = String::new();
    io::stdin().read_line(&mut input).unwrap();
    let nums: Vec<i64> = input.trim().split_whitespace().map(|x| x.parse().unwrap()).collect();
    let (n, r, g, b) = (nums[0], nums[1], nums[2], nums[3]);

    let mut pow2 = vec![1i64; 10_001];
    for i in 1..pow2.len() {
        pow2[i] = pow2[i - 1] * 2 % MOD;
    }

    let mut fact = vec![1i64; 2 * 100_000 + 1];
    for i in 1..fact.len() {
        fact[i] = fact[i - 1] * i as i64 % MOD;
    }

    let mut fact_inv = vec![0i64; fact.len()];
    fact_inv[fact.len() - 1] = mod_pow(fact[fact.len() - 1], MOD - 2, MOD);
    for i in (1..fact.len()).rev() {
        fact_inv[i - 1] = fact_inv[i] * i as i64 % MOD;
    }

    fn combi(a: i64, b: i64, fact: &Vec<i64>, fact_inv: &Vec<i64>) -> i64 {
        if b <= a {
            fact[a as usize] * fact_inv[b as usize] % MOD * fact_inv[(a - b) as usize] % MOD
        } else {
            0
        }
    }

    let mut ans: i64 = 0;

    for x in (0..=r.min(g)).rev() {
        for y in (0..=g.min(b)).rev() {
            for z in (0..=b.min(r)).rev() {
                if x * 2 + y * 2 + z * 2 > n {
                    break;
                }

                let r_rest = r - x - z;
                let g_rest = g - x - y;
                let b_rest = b - y - z;

                if r_rest < 0 || g_rest < 0 || b_rest < 0 {
                    continue;
                }

                let need = x + y + z + r_rest + g_rest + b_rest - 1;

                if x * 2 + y * 2 + z * 2 + r_rest + g_rest + b_rest + need > n {
                    break;
                }

                let empty = n - x * 2 - y * 2 - z * 2 - r_rest - g_rest - b_rest - need;

                let ax = combi(x + y + z + r_rest + g_rest + b_rest, x, &fact, &fact_inv)
                    * combi(y + z + r_rest + g_rest + b_rest, y, &fact, &fact_inv) % MOD
                    * combi(z + r_rest + g_rest + b_rest, z, &fact, &fact_inv) % MOD
                    * combi(r_rest + g_rest + b_rest, r_rest, &fact, &fact_inv) % MOD
                    * combi(g_rest + b_rest, g_rest, &fact, &fact_inv) % MOD;

                let ay = pow2[(x + y + z) as usize];

                let az = combi(empty + (x + y + z + r_rest + g_rest + b_rest), empty, &fact, &fact_inv);

                ans = (ans + ax * ay % MOD * az % MOD) % MOD;
            }
        }
    }

    println!("{}", ans);
}