/*
 * RiskScape™ Copyright New Zealand Institute for Earth Science Limited
 * (Earth Sciences New Zealand) is distributed for research purposes only
 * under the terms of AGPLv3.
 *
 * RiskScape™ Copyright 2026 New Zealand Institute for Earth Science
 * Limited (Earth Sciences New Zealand). All rights reserved. Source code
 * available under the AGPLv3.
 * 
 * This program is free software: you can redistribute it and/or modify it under
 *  the terms of the GNU Affero General Public License as published by the Free
 *  Software Foundation, either version 3 of the License, or (at your option) any
 *  later version.
 * 
 * This program is distributed for RESEARCH PURPOSES ONLY, in the hope that it will
 * be useful for research and education initiatives.
 * 
 * If you are not a researcher, or you are a researcher who wishes to use this
 * program on terms other than AGPLv3 (including those who wish to restrict the
 * distribution of any source code created using this program), please contact:
 * https://riskscape.org.nz
 * 
 * This program is distributed WITHOUT ANY WARRANTY; without even the implied
 * warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * Affero General Public License for more details.  You should have received a copy
 * of the GNU Affero General Public License along with this program.  If not, see
 * <http://www.gnu.org/licenses/>.
 * 
 * By way of summary only, under the AGPLv3:
 *     • Permissions of this strongest copyleft license are conditioned
 *       on making available complete source code of licensed works and
 *       modifications, which include larger works using a licensed work,
 *       under the same license.
 *     • Copyright and license notices must be preserved.
 *     • Contributors provide an express grant of patent rights.
 *     • When a modified version is used to provide a service over a
 *       network, the complete source code of the modified version must be made
 *       available.
 */
package nz.org.riskscape.engine.raster;

import static org.hamcrest.Matchers.*;
import static org.junit.Assert.*;

import java.util.Arrays;
import java.util.List;

import javax.media.jai.PlanarImage;

import org.geotools.coverage.grid.GridCoverage2D;
import org.geotools.geometry.Position2D;
import org.geotools.geometry.jts.ReferencedEnvelope;
import org.geotools.coverage.util.CoverageUtilities;
import org.hamcrest.Description;
import org.hamcrest.Matcher;
import org.hamcrest.TypeSafeMatcher;
import org.junit.Test;
import org.junit.runner.RunWith;
import org.junit.runners.Parameterized;
import org.junit.runners.Parameterized.Parameters;
import org.geotools.api.geometry.Position;

import org.locationtech.jts.geom.CoordinateSequence;
import org.locationtech.jts.geom.GeometryFactory;
import org.locationtech.jts.geom.LinearRing;
import org.locationtech.jts.geom.Polygon;
import org.locationtech.jts.geom.impl.CoordinateArraySequence;

import it.geosolutions.jaiext.range.NoDataContainer;
import nz.org.riskscape.engine.ProjectTest;
import nz.org.riskscape.engine.Tuple;
import nz.org.riskscape.engine.query.TupleUtils;
import nz.org.riskscape.engine.raster.VectorToRaster.DrawFeatureResult;
import nz.org.riskscape.engine.relation.ListRelation;
import nz.org.riskscape.engine.relation.SpatialMetadata;
import nz.org.riskscape.engine.rl.DefaultOperators;
import nz.org.riskscape.engine.types.Struct;
import nz.org.riskscape.engine.types.Types;
import nz.org.riskscape.rl.ast.PropertyAccess;

@RunWith(Parameterized.class)
public class VectorToRasterTest extends ProjectTest {

  @Parameters
  public static List<Double> scaleFactors() {
    return Arrays.asList(
        1.0D,
        2.0D,
        0.5D,
        0.75D
    );
  }

  Double scaleFactor;

  GeometryFactory gf = project.getSridSet().getGeometryFactory(project.getDefaultCrs());

  Struct type = Struct.of("value", Types.INTEGER, "geom", Types.GEOMETRY);

  ListRelation relation = new ListRelation(
      // simple poly
      Tuple.ofValues(type, 1L, gf.createPolygon(polybox(3, 3, 8, 8))),
      Tuple.ofValues(type, 11L, gf.createPolygon(polybox(1, 1, 2, 2))),

      // multi poly
      Tuple.ofValues(type, 2L,
        gf.createMultiPolygon(new Polygon[] {
          gf.createPolygon(polybox(10, 10, 20, 20)),
          gf.createPolygon(polybox(30, 10, 40, 20))
        })
      ),

      // poly with holes
      Tuple.ofValues(type, 3L,
        gf.createPolygon(
            gf.createLinearRing(polybox(50, 10, 200, 100)),
            new LinearRing[] {
                gf.createLinearRing(polybox(60, 70, 70, 80)),
                gf.createLinearRing(polybox(150, 50, 170, 60))
            }
        )
      ),

      // concentric multi polys with holes
      Tuple.ofValues(type, 4L,
          gf.createMultiPolygon(new Polygon[] {
              gf.createPolygon(
                  gf.createLinearRing(polybox(5, 50, 45, 190)),
                  new LinearRing[] {
                      gf.createLinearRing(polybox(10, 55, 40, 185)),
                  }
              ),
              gf.createPolygon(
                  gf.createLinearRing(polybox(15, 60, 35, 120)),
                  new LinearRing[] {
                      gf.createLinearRing(polybox(20, 65, 30, 90)),
                  }
              ),
          })
      )
  ).withSpatialMetadata(new SpatialMetadata(project.getDefaultCrs(), TupleUtils.findRequiredGeometryMember(type)));
  // test the different shapes, multi and single
  // test nzgd and nztm

  GridCoverage2D coverage;
  double boundsWidth = 200;
  double boundsHeight = 200;

  double boundsX = 0;
  double boundsY = 0;


  public VectorToRasterTest(Double scaleFactor) {
    this.scaleFactor = scaleFactor;

    this.project.getFunctionSet().insertFirst(DefaultOperators.INSTANCE);
  }


  @Test
  public void smokeTest() throws Exception {
    relation = new ListRelation(
        Tuple.ofValues(type, 100L, gf.createPolygon(polybox(1, 1, 4, 4))))
        .withSpatialMetadata(new SpatialMetadata(project.getDefaultCrs(), TupleUtils.findRequiredGeometryMember(type)));
    boundsHeight = 5;
    boundsWidth = 5;

    draw();

    if (scaleFactor < 1D) {
      assertThat(coverage, hasDataAtPoint(100, 0.5, 0.5));
    } else {
      assertThat(coverage, hasNoDataAtPoint(0.5, 0.5));
    }
    assertThat(coverage, hasDataAtPoint(100, 2, 2));
  }

  @Test
  public void checkNoDataValues() throws Exception {
    draw();

    // four corners
    if (scaleFactor >= 1D) {
      assertThat(coverage, hasNoDataAtPoint(0.5, 0.5));
    }
    assertThat(coverage, hasNoDataAtPoint(199.5, 199.5));
    assertThat(coverage, hasNoDataAtPoint(0.5, 199.5));
    assertThat(coverage, hasNoDataAtPoint(199.5, 0.5));

    // in between multi poly
    assertThat(coverage, hasNoDataAtPoint(25.5, 15.5));

    // in the holes of the poly
    assertThat(coverage, hasNoDataAtPoint(65.5, 75.5));
    assertThat(coverage, hasNoDataAtPoint(160.5, 55.5));

    // in the holes of the multi poly
    assertThat(coverage, hasNoDataAtPoint(20.5, 150));
    assertThat(coverage, hasNoDataAtPoint(25.5, 70.5));
  }

  @Test
  public void checkDataSampling() throws Exception {
    draw();

    // simple poly
    assertThat(coverage, hasDataAtPoint(1, 6.5, 6.5));
    // a 1x1 polygon
    assertThat(coverage, hasDataAtPoint(11, 1.5, 1.5));

    // multi poly
    assertThat(coverage, hasDataAtPoint(2, 15.5, 15.5));

    assertThat(coverage, hasDataAtPoint(3, 51.5, 11.5));


    assertThat(coverage, hasDataAtPoint(4, 7.5, 53.5));
    assertThat(coverage, hasDataAtPoint(4, 17.5, 62.5));
  }

  @Test
  public void featuresOutsideOfBoundsAreNotRendered() throws Exception {
    boundsX = 10;
    boundsY = 10;

    draw();

    // multi poly is in
    assertThat(coverage, hasDataAtPoint(2, 15.5, 15.5));

    // the pixel value of 1 should never occur
    PlanarImage image = (PlanarImage) coverage.getRenderedImage();
    int width = image.getWidth();
    int height = image.getHeight();
    float[] pixelHolder = new float[] {Float.NaN};
    for (int x = image.getMinX(); x < width; x++) {
      for (int y = image.getMinY(); y < height; y++) {
          image.getTile(image.XToTileX(x), image.YToTileY(y)).getPixel(x, y, pixelHolder);
          assertNotEquals(pixelHolder[0], 1.0F, 0.0001);
      }
    }
  }

  @Test
  public void featuresOnTheBorderOfBoundsAreRendered() throws Exception {
    boundsX = 4;
    boundsY = 4;
    boundsWidth = 3;
    boundsHeight = 3;

    draw();

    // simple poly is in
    assertThat(coverage, hasDataAtPoint(1, 6.5, 6.5));
  }

  @Test
  public void pointsOnTheEdgesOfTheBoundsAreCovered() throws Exception {
    relation = new ListRelation(
        Tuple.ofValues(type, 1L, gf.createPoint(coords(1, 1))),
        Tuple.ofValues(type, 2L, gf.createPoint(coords(6, 6)))
    ).inferSpatialMetadata(project);

    draw(relation.calculateBounds().get());

    // due to rounding during pixel transformation, the exact points of the bounds are not reliably translated
    // so we add a small fraction on make sure we stay within the bounds during crs -> grid transformation
    hasDataAtPoint(1, 1.0001, 1.0001);
    hasDataAtPoint(2, 5.9999, 5.9999);
  }

  @Test
  public void drawsFeaturesFromTuples() throws Exception {
    ReferencedEnvelope bounds = new ReferencedEnvelope(
        boundsX, boundsX + boundsWidth, boundsY, boundsY + boundsHeight,
        project.getDefaultCrs()
    );
    VectorToRaster v2r = new VectorToRaster();
    v2r.initialize(bounds, scaleFactor, VectorToRaster.PixelStrategy.OVERWRITE);

    // null geometry
    assertThat(v2r.drawFeature(2F, null), is(DrawFeatureResult.SKIPPED_NO_VALUE_OR_GEOMETRY));
    // null value
    assertThat(
        v2r.drawFeature((Number)null, gf.createPolygon(polybox(0, 0, 1, 1))),
        is(DrawFeatureResult.SKIPPED_NO_VALUE_OR_GEOMETRY)
    );
    assertThat(
        v2r.drawFeature(11D, gf.createPolygon(polybox(1, 1, 2, 2))),
        is(DrawFeatureResult.DRAWN)
    );
    assertThat(
        v2r.drawFeature(1L, gf.createPolygon(polybox(3, 3, 8, 8))),
        is(DrawFeatureResult.DRAWN)
    );

    // out of bounds
    assertThat(
        v2r.drawFeature(1L, gf.createPolygon(polybox(-100, -100, -80, -80))),
        is(DrawFeatureResult.OUT_OF_BOUNDS)
    );

    coverage = v2r.constructCoverage("test");
    if (scaleFactor < 1D) {
      // at scale factors below 1 polybox(1, 1, 2, 2) moves into this pixel
      assertThat(coverage, hasDataAtPoint(11, 0.5, 0.5));
    } else {
      assertThat(coverage, hasNoDataAtPoint(0.5, 0.5));
    }

    // a 1x1 polygon
    assertThat(coverage, hasDataAtPoint(11, 1.5, 1.5));

    assertThat(coverage, hasDataAtPoint(1, 6.5, 6.5));
  }

  private Matcher<GridCoverage2D> hasDataAtPoint(double expectedValue, double x, double y) {
    return new TypeSafeMatcher<GridCoverage2D>() {

      @Override
      public void describeTo(Description description) {
        description.appendText(String.format("%f at (%f %f)", expectedValue, x, y));
      }

      @Override
      protected boolean matchesSafely(GridCoverage2D item) {
        Position2D point = new Position2D(project.getDefaultCrs(), x, y);

        double[] container = new double[1];
        item.evaluate((Position) point, container);

        return container[0] == expectedValue;
      }

      @Override
      protected void describeMismatchSafely(GridCoverage2D item, Description mismatchDescription) {
        NoDataContainer noData = CoverageUtilities.getNoDataProperty(coverage);
        Position2D point = new Position2D(project.getDefaultCrs(), x, y);
        double noDataValue = noData.getAsSingleValue();
        double[] container = new double[1];
        item.evaluate((Position) point, container);
        String value = container[0] == noDataValue ? "NODATA" : Double.toString(container[0]);

        mismatchDescription.appendText(String.format("got %s", value));
      }
    };
  }

  private void draw() throws Exception {
    ReferencedEnvelope bounds = new ReferencedEnvelope(boundsX, boundsX + boundsWidth, boundsY, boundsY + boundsHeight,
        project.getDefaultCrs());
    draw(bounds);
  }

  private void draw(ReferencedEnvelope bounds) throws Exception {
    VectorToRaster v2r = new VectorToRaster();

    coverage = v2r.convert(
        relation,
        expressionRealizer.realize(relation.getType(), PropertyAccess.of("value")).get(),
        scaleFactor,
        bounds,
        "test");
  }

  private Matcher<GridCoverage2D> hasNoDataAtPoint(double x, double y) {
    return new TypeSafeMatcher<GridCoverage2D>() {

      @Override
      public void describeTo(Description description) {
        description.appendText(String.format("no data at %f %f", x, y));
      }

      @Override
      protected boolean matchesSafely(GridCoverage2D item) {
        NoDataContainer noData = CoverageUtilities.getNoDataProperty(coverage);
        Position2D point = new Position2D(project.getDefaultCrs(), x, y);

        double[] container = new double[1];
        item.evaluate((Position) point, container);


        return noData.getAsRange().contains(container[0]);
      }

      @Override
      protected void describeMismatchSafely(GridCoverage2D item, Description mismatchDescription) {
        NoDataContainer noData = CoverageUtilities.getNoDataProperty(coverage);
        Position2D point = new Position2D(project.getDefaultCrs(), x, y);
        double noDataValue = noData.getAsSingleValue();
        double[] container = new double[1];
        item.evaluate((Position) point, container);
        String value = container[0] == noDataValue ? "NODATA" : Double.toString(container[0]);

        mismatchDescription.appendText(String.format("got %s", value));
      }
    };
  }

  private CoordinateSequence polybox(double x1, double y1, double x2, double y2) {
    return coords(x1, y1, x2, y1, x2, y2, x1, y2, x1, y1);
  }

  private CoordinateSequence coords(double... xThenYs) {
    if ((xThenYs.length % 2) != 0) {
      throw new IllegalArgumentException("uneven number of points");
    }

    CoordinateArraySequence sequence = new CoordinateArraySequence(xThenYs.length / 2);
    for (int arrayIdx = 0, seqIdx = 0; arrayIdx < xThenYs.length; seqIdx++) {
      sequence.setOrdinate(seqIdx, 0, xThenYs[arrayIdx++]);
      sequence.setOrdinate(seqIdx, 1, xThenYs[arrayIdx++]);
    }

    return sequence;
  }

}
