#pragma GCC target ("avx2")
#pragma GCC optimize ("unroll-loops")
#pragma GCC optimize ("O3")
#include "bits/stdc++.h"
#include <unordered_set>
#include <unordered_map>
#include <random>
using namespace std;
typedef long long ll;
const ll MOD = /*1'000'000'007LL;*/ 998'244'353LL;
#define pb push_back
#define mp make_pair
#define all(x) (x).begin(), (x).end()
#define rep(i, n) for(int (i)=0; (i)<(n); (i)++)
const int dx[4]={ 1,0,-1,0 };
const int dy[4]={ 0,1,0,-1 };

template<long long Mod>
struct mint{
    long long x;
    mint(long long x=0) :x((x% Mod+Mod)%Mod){}

    mint operator-() const { return mint(-x); }

    mint& operator+=(const mint a) {
        if((x += a.x) >= Mod) x -= Mod;
        return *this;
    }
    mint& operator-=(const mint a) {
        if((x += Mod-a.x) >= Mod) x -= Mod;
        return *this;
    }
    mint& operator*=(const mint a) {
        (x *= a.x) %= Mod;
        return *this;
    }

    mint operator^(long long t) const {
        mint ret(1);
        mint a = *this;
        while(t){
            if(t&1){
                ret *= a;
            }
            a *= a;
            t >>= 1;
        }
        return ret;
    }
    mint inv() const { return *this ^ (Mod-2); }
    mint& operator/=(const mint a) { return *this *= a.inv(); }

    mint operator+(const mint a) const { return mint(*this) += a; }
    mint operator-(const mint a) const { return mint(*this) -= a; }
    mint operator*(const mint a) const { return mint(*this) *= a; }
    mint operator/(const mint a) const { return mint(*this) /= a; }

    friend istream& operator>>(istream& is, mint& a) { return is >> a.x; }
    friend ostream& operator<<(ostream& os, const mint& a) { return os << a.x; }
};

template<typename T>
struct Combination{
private:
    int N;
    vector<T> f, i;
public:
    void init(int n){
        N = n;
        f.resize(N+1); f[0] = T(1);
        for(int j=1; j<=N; j++){
            f[j] = f[j-1] * T(j);
        }
        i.resize(N+1); i[N] = T(1) / f[N];
        for(int j=N-1; j>=0; j--){
            i[j] = i[j+1] * T(j+1);
        }
    }
    T fact(int n){
        return f[n];
    }
    T P(int n, int r){
        return f[n] * i[n-r];
    }
    T C(int n, int r){
        return f[n] * i[r] * i[n-r];
    }
};

const long long conv_mod = 998'244'353LL;
const mint<conv_mod> g = 3;

template<typename T>
void DFT(vector<T>& f){
    if(f.size() == 1) return;
    int n = f.size();
    vector<T> f0(n>>1), f1(n>>1);
    for(int i=0; i<n; i++){
        if((i&1) == 0) f0[i>>1] = f[i];
        else f1[i>>1] = f[i];
    }
    DFT(f0);
    DFT(f1);
    T zeta = g ^ ((conv_mod-1) / f.size());
    T zeta_i = 1;
    for(int i=0; i<n; i++){
        f[i] = f0[i % (n/2)] + zeta_i*f1[i % (n/2)];
        zeta_i *= zeta;
    }
}

template<typename T>
void invDFT(vector<T>& f){
    if(f.size() == 1) return;
    int n = f.size();
    vector<T> f0(n>>1), f1(n>>1);
    for(int i=0; i<n; i++){
        if((i&1) == 0) f0[i>>1] = f[i];
        else f1[i>>1] = f[i];
    }
    invDFT(f0);
    invDFT(f1);
    T inv_zeta = T(1) / (g ^ ((conv_mod-1) / f.size()));
    T zeta_i = 1;
    for(int i=0; i<n; i++){
        f[i] = f0[i % (n/2)] + zeta_i*f1[i % (n/2)];
        zeta_i *= inv_zeta;
    }
}

vector<mint<conv_mod>> convolution(vector<mint<conv_mod>> A, vector<mint<conv_mod>> B){
    int n = 1;
    while(n < A.size() + B.size() - 1) n <<= 1;
    vector<mint<conv_mod>> a(n), b(n), c(n);
    for(int i=0; i<A.size(); i++) a[i] = A[i];
    for(int i=0; i<B.size(); i++) b[i] = B[i];
    DFT(a);
    DFT(b);
    for(int i=0; i<n; i++) c[i] = a[i] * b[i];
    invDFT(c);
    vector<mint<conv_mod>> ret(A.size() + B.size() - 1);
    for(int i=0; i<ret.size(); i++) ret[i] = c[i] / mint<conv_mod>(n);
    return ret;
}

int N;
int h[100000];
int c[100000];

Combination<mint<MOD>> comb;
mint<MOD> memo[100001];

vector<mint<MOD>> solve(int l, int r){
    if(r-l == 1){
        vector<mint<MOD>> ret(c[l]/2+1);
        rep(i, c[l]/2+1){
            ret[i] = comb.C(c[l], 2*i) * memo[i];
        }
        return ret;
    }
    int m = (l+r)/2;
    auto L = solve(l, m);
    auto R = solve(m, r);
    return convolution(L, R);
}

signed main(){
    cin >> N;
    rep(i, 2*N){
        cin >> h[i];
        h[i]--;
        c[h[i]]++;
    }

    comb.init(2*N);
    memo[0] = 1;
    rep(i, N) memo[i+1] = memo[i] * (2*i+1);
    auto v = solve(0, 100000);

    mint<MOD> ans = 0;
    rep(i, v.size()){
        if(i%2 == 0) ans += memo[N-i] * v[i];
        else ans -= memo[N-i] * v[i];
    }
    cout << ans << endl;
}

Measured memory_bytes: 9944