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

using ll = long long;
using ld = long double;

#define all(a) begin(a), end(a)
#define len(a) int((a).size())

#ifdef LOCAL
    #include "debug.h"
#else
    #define dbg(...)
    #define dprint(...)
    #define debug if constexpr (false)
    #define draw_tree(...)
    #define draw_array(...)
#endif // LOCAL

/*
 ! WARNING: MOD must be prime.
 ! WARNING: MOD must be less than 2^30.
 * Use .get() to get the stored value.
 */
template<uint32_t mod>
class montgomery {
    static_assert(mod < uint32_t(1) << 30, "mod < 2^30");
    using mint = montgomery<mod>;

private:
    uint32_t value;

    static constexpr uint32_t inv_neg_mod = []() {
        uint32_t x = mod;
        for (int i = 0; i < 4; ++i) {
            x *= uint32_t(2) - mod * x;
        }
        return -x;
    }();
    static_assert(mod * inv_neg_mod == -1);

    static constexpr uint32_t neg_mod = (-uint64_t(mod)) % mod;

    static uint32_t reduce(const uint64_t &value) {
        return (value + uint64_t(uint32_t(value) * inv_neg_mod) * mod) >> 32;
    }

    inline static const mint ONE = mint(1);

public:
    montgomery() : value(0) {}
    montgomery(const mint &x) : value(x.value) {}

    template<typename T, typename U = std::enable_if_t<std::is_integral<T>::value>>
    montgomery(const T &x) : value(!x ? 0 : reduce(int64_t(x % int32_t(mod) + int32_t(mod)) * neg_mod)) {}

    static constexpr uint32_t get_mod() {
        return mod;
    }

    uint32_t get() const {
        auto real_value = reduce(value);
        return real_value < mod ? real_value : real_value - mod;
    }

    template<typename T>
    mint power(T degree) const {
        degree = (degree % int32_t(mod - 1) + int32_t(mod - 1)) % int32_t(mod - 1);
        mint prod = 1, a = *this;
        for (; degree > 0; degree >>= 1, a *= a)
            if (degree & 1)
                prod *= a;

        return prod;
    }

    mint inv() const {
        return power(-1);
    }

    mint& operator=(const mint &x) {
        value = x.value;
        return *this;
    }

    mint& operator+=(const mint &x) {
        if (int32_t(value += x.value - (mod << 1)) < 0) {
            value += mod << 1;
        }
        return *this;
    }

    mint& operator-=(const mint &x) {
        if (int32_t(value -= x.value) < 0) {
            value += mod << 1;
        }
        return *this;
    }

    mint& operator*=(const mint &x) {
        value = reduce(uint64_t(value) * x.value);
        return *this;
    }

    mint& operator/=(const mint &x) {
        return *this *= x.inv();
    }

    friend mint operator+(const mint &x, const mint &y) {
        return mint(x) += y;
    }

    friend mint operator-(const mint &x, const mint &y) {
        return mint(x) -= y;
    }

    friend mint operator*(const mint &x, const mint &y) {
        return mint(x) *= y;
    }

    friend mint operator/(const mint &x, const mint &y) {
        return mint(x) /= y;
    }

    mint& operator++() {
        return *this += ONE;
    }

    mint& operator--() {
        return *this -= ONE;
    }

    mint operator++(int) {
        mint prev = *this;
        *this += ONE;
        return prev;
    }

    mint operator--(int) {
        mint prev = *this;
        *this -= ONE;
        return prev;
    }

    mint operator-() const {
        return mint(0) - *this;
    }

    bool operator==(const mint &x) const {
        return get() == x.get();
    }

    bool operator!=(const mint &x) const {
        return get() != x.get();
    }

    bool operator<(const mint &x) const {
        return get() < x.get();
    }

    template<typename T>
    explicit operator T() {
        return get();
    }

    friend std::istream& operator>>(std::istream &in, mint &x) {
        std::string s;
        in >> s;
        x = 0;
        bool neg = s[0] == '-';
        for (const auto c : s)
            if (c != '-')
                x = x * 10 + (c - '0');

        if (neg)
            x *= -1;

        return in;
    }

    friend std::ostream& operator<<(std::ostream &out, const mint &x) {
        return out << x.get();
    }

    static int32_t primitive_root() {
        if constexpr (mod == 1'000'000'007)
            return 5;
        if constexpr (mod == 998'244'353)
            return 3;
        if constexpr (mod == 786433)
            return 10;

        static int root = -1;
        if (root != -1)
            return root;

        std::vector<int> primes;
        int value = mod - 1;
        for (int i = 2; i * i <= value; i++)
            if (value % i == 0) {
                primes.push_back(i);
                while (value % i == 0)
                    value /= i;
            }

        if (value != 1)
            primes.push_back(value);

        for (int r = 2;; r++) {
            bool ok = true;
            for (auto p : primes)
                if ((mint(r).power((mod - 1) / p)).get() == 1) {
                    ok = false;
                    break;
                }

            if (ok)
                return root = r;
        }
    }
};

constexpr uint32_t MOD = 1'000'000'007;
// constexpr uint32_t MOD = 998'244'353;
using mint = montgomery<MOD>;

/*
 ! WARNING: MOD must be prime.
 * Define modular int class above it.
 * No need to run any init function, it dynamically resizes the data.
 */
namespace combinatorics {
    std::vector<mint> fact_, ifact_, inv_;

    void resize_data(int size) {
        if (fact_.empty()) {
            fact_ = {mint(1), mint(1)};
            ifact_ = {mint(1), mint(1)};
            inv_ = {mint(0), mint(1)};
        }
        for (int pos = fact_.size(); pos <= size; pos++) {
            fact_.push_back(fact_.back() * mint(pos));
            inv_.push_back(-inv_[MOD % pos] * mint(MOD / pos));
            ifact_.push_back(ifact_.back() * inv_[pos]);
        }
    }

    struct combinatorics_info {
        std::vector<mint> &data;

        combinatorics_info(std::vector<mint> &data) : data(data) {}

        mint operator[](int pos) {
            if (pos >= static_cast<int>(data.size())) {
                resize_data(pos);
            }
            return data[pos];
        }
    } fact(fact_), ifact(ifact_), inv(inv_);

    // From n choose k.
    // O(max(n)) in total.
    mint choose(int n, int k) {
        if (n < k || k < 0 || n < 0) {
            return mint(0);
        }
        return fact[n] * ifact[k] * ifact[n - k];
    }

    // From n choose k.
    // O(min(k, n - k)).
    mint choose_slow(int64_t n, int64_t k) {
        if (n < k || k < 0 || n < 0) {
            return mint(0);
        }
        k = std::min(k, n - k);
        mint result = 1;
        for (int i = k; i >= 1; i--) {
            result *= (n - i + 1);
            result *= inv[i];
        }
        return result;
    }

    // Number of balanced bracket sequences with n open and m closing brackets.
    mint catalan(int n, int m) {
        if (m > n || m < 0) {
            return mint(0);
        }
        return choose(n + m, m) - choose(n + m, m - 1);
    }

    // Number of balanced bracket sequences with n open and closing brackets.
    mint catalan(int n) {
        return catalan(n, n);
    }
} // namespace combinatorics

using namespace combinatorics;

vector<mint> Lagrange(vector<mint> &x, vector<mint> &y) {
  int n = x.size();
  vector<mint> ans(n, 0);
  #undef all
  vector<mint> all(n + 1, 0); //(x - x0) * (x - x1) * ... * (x - x(n-1))
  all[0] = 1;
  for (int i = 0; i < n; i++) {
    for (int j = n; j >= 0; j--) {
      all[j] *= -x[i];
      if (j) all[j] += all[j - 1];
    }
  }
  for (int i = 0; i < n; i++) {
    vector <mint> up(n); //all / (x - xi)
    mint rem = all[n];
    for (int j = n - 1; j >= 0; --j) {
      up[j] = rem;
      rem = all[j] + rem * x[i];
    }
    mint down = 1;
    for (int j = 0; j < n; j++) if (i != j) down *= x[i] - x[j];
    up.resize(n); down = down.inv() * y[i];
    for (int j = 0; j < n; j++) ans[j] += up[j] * down;
  }
  return ans;
}

int main() {
    cin.tie(nullptr)->sync_with_stdio(false);

    int A, B, C, x;
    cin >> A >> B >> C >> x;

    const int S = A + B + C;
    if (S % 2 != 0) {
        cout << "0\n";
        return 0;
    }

    vector<vector<mint>> dp(S + 1);

    vector<mint> cat(S + 1);
    for (int n = 0; n <= S; n += 2) {
        cat[n] = catalan(n / 2);
    }

    dp[0] = {1};
    for (int n = 2; n <= S; n += 2) {
        dp[n].resize(n / 2 + 1);
        for (int prev = 0; prev < n; prev += 2) {
            for (int cnt = 1; cnt <= prev / 2 + 1; cnt++) {
                dp[n][cnt] += dp[prev][cnt - 1] * cat[n - prev - 2];
            }
        }
    }

    vector<vector<mint>> vals(S + 1, vector<mint>(S / 2 + 2));
    for (int n = 0; n <= S; n += 2) {
        for (int cnt = 0; cnt <= n / 2; cnt++) {
            dp[n][cnt] *= ifact[cnt];
        }
        for (int x = 0; x < len(vals[n]); x++) {
            for (int i = n / 2; i >= 0; i--) {
                (vals[n][x] *= x) += dp[n][i];
            }
        }
    }

    vector<mint> tot(S / 2 + 2);
    for (int a = 0; a <= A; a++) {
        for (int b = 0; b <= B; b++) {
            int c = S / 2 - a - b;
            if (c < 0) {
                continue;
            }
            if (c > C) {
                continue;
            }

            mint coeff = 0;
            for (int cb = 0; cb <= a && cb <= B - b; cb++) {
                coeff += choose(a, cb) * choose(c, B - b - cb) * choose(b, A - a - (c - (B - b - cb)));
            }
            for (int i = 0; i < len(tot); i++) {
                tot[i] += coeff * vals[2 * a][i] * vals[2 * b][i] * vals[2 * c][i];
            }
        }
    }

    vector<mint> xs(len(tot));
    iota(xs.begin(), xs.end(), 0);
    tot = Lagrange(xs, tot);

    mint ans = 0;
    for (int i = 1; i < len(tot); i++) {
        ans += tot[i] * mint(i).power(x) * fact[i];
    }
    cout << ans << '\n';
}