package cz.vutbr.feec.utko.xonder02.lda;

import java.util.Random;

public class Compute {

    public static final double EULER = 0.577215664901532860606512090082402431;

    static double psiCoeficients[] = {
        8.33333333333333333333E-2,
        -2.10927960927960927961E-2,
        7.57575757575757575758E-3,
        -4.16666666666666666667E-3,
        3.96825396825396825397E-3,
        -8.33333333333333333333E-3,
        8.33333333333333333333E-2
    };

    public static double polevl(double x, double[] coef) {
        double ans = coef[0];

        for (int i = 1; i < coef.length; i++) {
            ans = ans * x + coef[i];
        }

        return ans;
    }

    public static double psi(double x) {
        double p, q, nz, s, w, y, z;
        int i, n;
        boolean negative = false;

        nz = 0.0;

        if (x <= 0.0) {
            negative = true;
            q = x;
            p = Math.floor(q);
            if (p == q) {
                return Double.POSITIVE_INFINITY;
            }
            /* Remove the zeros of tan(NPY_PI x)
             * by subtracting the nearest integer from x
             */
            nz = q - p;
            if (nz != 0.5) {
                if (nz > 0.5) {
                    p += 1.0;
                    nz = q - p;
                }
                nz = Math.PI / Math.tan(Math.PI * nz);
            }
            else {
                nz = 0.0;
            }
            x = 1.0 - x;
        }

        /* check for positive integer up to 10 */
       if ((x <= 10.0) && (x == Math.floor(x))) {
            y = 0.0;
            n = (int)x;
            for (i = 1; i < n; i++) {
                w = i;
                y += 1.0 / w;
            }
            y -= EULER;

            if (negative) {
                y -= nz;
            }

            return y;
        }

        s = x;
        w = 0.0;
        while (s < 10.0) {
            w += 1.0 / s;
            s += 1.0;
        }

        if (s < 1.0e17) {
            z = 1.0 / (s * s);
            y = z * polevl(z, psiCoeficients);
        } else {
            y = 0.0;
        }

        y = Math.log(s) - (0.5 / s) - y - w;

        if (negative) {
            y -= nz;
        }

        return y;
    }

    public static double[][] psi(double[][] input) {
        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = psi(input[i][j]);
            }
        }

        return output;
    }

    public static double[][] dirichletExpectation(double[][] input) {
        return subPerRow(psi(input), psi(new double[][] { sumPerRow(input) })[0]);
    }

    public static double sum(double[][] input) {
        double sum = .0;

        for (double[] anInput : input) {
            for (int j = 0; j < input[0].length; j++) {
                sum += anInput[j];
            }
        }

        return sum;
    }

    public static double sum(double[] input) {
        double sum = .0;

        for (double anInput : input) {
            sum += anInput;
        }

        return sum;
    }

    public static double[] sumPerRow(double[][] input) {
        double[] output = new double[input.length];

        for (int i = 0; i < input.length; i++) {
            double sum = 0;
            for (int j = 0; j < input[0].length; j++) {
                sum += input[i][j];
            }
            output[i] = sum;
        }

        return output;
    }

    public static double[][] add(double[][] input, double add) {
        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = input[i][j] + add;
            }
        }

        return output;
    }

    public static double[][] sub(double[][] input, double substract) {
        return add(input, -substract);
    }

    public static double[][] sub(double[][] input, double[][] substract) throws Exception {
        checkSizes(input, substract);

        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = input[i][j] - substract[i][j];
            }
        }

        return output;
    }

    public static double[][] subAbsolute(double[][] input, double[][] substract) throws Exception {
        checkSizes(input, substract);

        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = Math.abs(input[i][j] - substract[i][j]);
            }
        }

        return output;
    }

    public static double[] subAbsolute(double[] input, double[] substract) {
        double[] output = new double[input.length];

        for (int i = 0; i < output.length; i++) {
            output[i] = Math.abs(input[i] - substract[i]);
        }

        return output;
    }

    public static double[][] subPerRow(double[][] input, double[] substract) {
        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = input[i][j] - substract[i];
            }
        }

        return output;
    }

    public static double[][] mul(double[][] input, double mul) {
        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = input[i][j] * mul;
            }
        }

        return output;
    }

    public static double[][] div(double[][] input, double[][] div) throws Exception {
        checkSizes(input, div);

        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = input[i][j] / div[i][j];
            }
        }

        return output;
    }

    public static double[][] div(double[][] input, double div) {
        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = input[i][j] / div;
            }
        }

        return output;
    }

    public static double[][] log(double[][] input) {
        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = Math.log(input[i][j]);
            }
        }

        return output;
    }

    public static double[][] exp(double[][] input) {
        double[][] output = new double[input.length][input[0].length];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = Math.exp(input[i][j]);
            }
        }

        return output;
    }

    public static double[][] expMutable(double[][] input) {
        for (int i = 0; i < input.length; i++) {
            for (int j = 0; j < input[0].length; j++) {
                input[i][j] = Math.exp(input[i][j]);
            }
        }

        return input;
    }

    public static double[][] addMutable(double[][] input, double[][] add) throws Exception {
        checkSizes(input, add);

        for (int i = 0; i < input.length; i++) {
            for (int j = 0; j < input[0].length; j++) {
                input[i][j] += add[i][j];
            }
        }

        return input;
    }

    public static double[][] addMutable(double[][] input, double add) {
        for (int i = 0; i < input.length; i++) {
            for (int j = 0; j < input[0].length; j++) {
                input[i][j] += add;
            }
        }

        return input;
    }

    public static double[][] addMutableForColumns(double[][] input, double[][] add, int[] columns) {
        for (int i = 0; i < input.length; i++) {
            int addColumn = 0;
            for (int column : columns) {
                input[i][column] += add[i][addColumn++];
            }
        }

        return input;
    }

    public static double[][] subMutable(double[][] input, double[][] sub) throws Exception {
        checkSizes(input, sub);

        for (int i = 0; i < input.length; i++) {
            for (int j = 0; j < input[0].length; j++) {
                input[i][j] -= sub[i][j];
            }
        }

        return input;
    }

    public static double[][] mulMutable(double[][] input, double mul) {
        for (int i = 0; i < input.length; i++) {
            for (int j = 0; j < input[0].length; j++) {
                input[i][j] *= mul;
            }
        }

        return input;
    }

    public static double[][] mulMutable(double[][] input, double[][] mul) throws Exception {
        checkSizes(input, mul);

        for (int i = 0; i < input.length; i++) {
            for (int j = 0; j < input[0].length; j++) {
                input[i][j] *= mul[i][j];
            }
        }

        return input;
    }

    public static double[][] absMutable(double[][] input) {
        for (int i = 0; i < input.length; i++) {
            for (int j = 0; j < input[0].length; j++) {
                input[i][j] = Math.abs(input[i][j]);
            }
        }

        return input;
    }

    public static double mean(double[][] input) {
        double sum = sum(input);
        return sum / (input.length * input[0].length);
    }

    public static double mean(double[] input) {
        double sum = sum(input);
        return sum / input.length;
    }

    public static double meanDiff(double[] x, double[] y) {
        double sum = sum(subAbsolute(x, y));
        return sum / x.length;
    }

    public static double[][] transpose(double[][] input) {
        double[][] output = new double[input[0].length][input.length];

        for (int i = 0; i < input.length; i++) {
            for (int j = 0; j < input[0].length; j++) {
                output[j][i] = input[i][j];
            }
        }

        return output;
    }

    public static double[][] outer(double[] x, double[] y) {
        double[][] output = new double[x.length][y.length];

        for (int i = 0; i < x.length; i++) {
            for (int j = 0; j < y.length; j++) {
                output[i][j] = x[i] * y[j];
            }
        }

        return output;
    }

    public static double[][] matrixMultiplication(double[][] x, double[][] y) {
        if (x[0].length != y.length) {
            throw new IllegalArgumentException();
        }

        double[][] output = new double[x.length][y[0].length];

        for (int i = 0; i < x.length; i++) {
            for (int j = 0; j < y[0].length; j++) {
                for (int k = 0; k < x[0].length; k++) {
                    output[i][j] += x[i][k] * y[k][j];
                }
            }
        }

        return output;
    }

    public static double[][] gammaDistribution(double shape, double scale, int x, int y) {
        double[][] output = new double[x][y];

        GammaDistribution gammaDistribution = new GammaDistribution(new Random(), shape, scale);

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = gammaDistribution.sample();
            }
        }

        return output;
    }

    public static double[][] ones(int x, int y) {
        double[][] output = new double[x][y];

        for (int i = 0; i < output.length; i++) {
            for (int j = 0; j < output[0].length; j++) {
                output[i][j] = 1.0;
            }
        }

        return output;
    }

    public static double[][] getColumns(double[][] matrix, int[] ids) {
        double[][] output = new double[matrix.length][ids.length];

        for (int i = 0; i < matrix.length; i++) {
            for (int j = 0; j < ids.length; j++) {
                output[i][j] = matrix[i][ids[j]];
            }
        }

        return output;
    }

    public static double[][] copy(double[][] input) {
        double [][] output = new double[input.length][input[0].length];

        for (int i = 0; i < input.length; i++) {
            System.arraycopy(input[i], 0, output[i], 0, input[0].length);
        }

        return output;
    }

    public static double cosineDistance(double[] a, double[] b) throws Exception {
        if (a.length != b.length) {
            throw new Exception("das");
        }

        double numerator = 0;
        double aSum = 0;
        double bSum = 0;

        for (int i = 0; i < a.length; i++) {
            numerator += a[i] * b[i];
            aSum += a[i] * a[i];
            bSum += b[i] * b[i];
        }

        return numerator / (Math.sqrt(aSum) * Math.sqrt(bSum));
    }

    private static void checkSizes(double[][] x, double[][] y) throws Exception {
        if (x.length != y.length || x[0].length != y[0].length) {
            throw new Exception("das");
        }
    }
}
