package sparseviewer.interp;

public class IdwInterpolator {

    public static InterpolationResult interpolate(double[] x,
                                                  double[] y,
                                                  double[] z,
                                                  int width,
                                                  int height,
                                                  double power) {
        Grid2D grid = createGrid(x, y, width, height);
        float[] out = new float[width * height];

        double exponent = power * 0.5;

        for (int iy = 0; iy < height; iy++) {
            double gy = grid.yAt(iy);

            for (int ix = 0; ix < width; ix++) {
                double gx = grid.xAt(ix);

                double numerator = 0.0;
                double denominator = 0.0;
                boolean exact = false;
                double exactValue = 0.0;

                for (int i = 0; i < x.length; i++) {
                    if (!isFinite(x[i]) || !isFinite(y[i]) || !isFinite(z[i])) {
                        continue;
                    }

                    double dx = x[i] - gx;
                    double dy = y[i] - gy;
                    double d2 = dx * dx + dy * dy;

                    if (d2 < 1e-24) {
                        exact = true;
                        exactValue = z[i];
                        break;
                    }

                    double w = 1.0 / Math.pow(d2, exponent);

                    numerator += w * z[i];
                    denominator += w;
                }

                if (exact) {
                    out[iy * width + ix] = (float) exactValue;
                } else if (denominator > 0.0) {
                    out[iy * width + ix] = (float) (numerator / denominator);
                } else {
                    out[iy * width + ix] = Float.NaN;
                }
            }
        }

        return new InterpolationResult(grid, out);
    }

    private static Grid2D createGrid(double[] x, double[] y, int width, int height) {
        double xmin = Double.POSITIVE_INFINITY;
        double xmax = Double.NEGATIVE_INFINITY;
        double ymin = Double.POSITIVE_INFINITY;
        double ymax = Double.NEGATIVE_INFINITY;

        for (int i = 0; i < x.length; i++) {
            if (isFinite(x[i])) {
                if (x[i] < xmin) xmin = x[i];
                if (x[i] > xmax) xmax = x[i];
            }

            if (isFinite(y[i])) {
                if (y[i] < ymin) ymin = y[i];
                if (y[i] > ymax) ymax = y[i];
            }
        }

        if (xmin == xmax) {
            xmin -= 1.0;
            xmax += 1.0;
        }

        if (ymin == ymax) {
            ymin -= 1.0;
            ymax += 1.0;
        }

        double xpad = 0.03 * (xmax - xmin);
        double ypad = 0.03 * (ymax - ymin);

        return new Grid2D(width, height,
                xmin - xpad,
                xmax + xpad,
                ymin - ypad,
                ymax + ypad);
    }

    private static boolean isFinite(double v) {
        return !Double.isNaN(v) && !Double.isInfinite(v);
    }
}