/*
 * Licensed to the Apache Software Foundation (ASF) under one or more
 * contributor license agreements.  See the NOTICE file distributed with
 * this work for additional information regarding copyright ownership.
 * The ASF licenses this file to You under the Apache License, Version 2.0
 * (the "License"); you may not use this file except in compliance with
 * the License.  You may obtain a copy of the License at
 *
 *     http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */
package org.apache.sis.referencing.operation.provider;

import java.awt.geom.AffineTransform;
import jakarta.xml.bind.annotation.XmlTransient;
import org.opengis.util.FactoryException;
import org.opengis.parameter.ParameterDescriptor;
import org.opengis.parameter.ParameterDescriptorGroup;
import org.opengis.parameter.ParameterValueGroup;
import org.opengis.referencing.cs.CartesianCS;
import org.opengis.referencing.cs.EllipsoidalCS;
import org.opengis.referencing.operation.CylindricalProjection;
import org.opengis.referencing.operation.MathTransform;
import org.opengis.referencing.operation.MathTransformFactory;
import org.apache.sis.parameter.Parameters;
import org.apache.sis.parameter.ParameterBuilder;
import org.apache.sis.metadata.iso.citation.Citations;
import org.apache.sis.referencing.operation.matrix.MatrixSIS;
import org.apache.sis.referencing.operation.transform.ContextualParameters;
import org.apache.sis.referencing.operation.transform.MathTransforms;
import org.apache.sis.referencing.util.j2d.ParameterizedAffine;
import org.apache.sis.referencing.internal.Resources;
import org.apache.sis.util.internal.DoubleDouble;
import org.apache.sis.util.internal.Constants;

import static java.lang.Math.*;


/**
 * The provider for <cite>"Equidistant Cylindrical (Spherical)"</cite> projection
 * (EPSG:1029, <span class="deprecated">EPSG:9823</span>).
 * In the particular case where the longitude of origin and the standard parallel are 0°,
 * this projection is also known as <cite>"Plate Carrée"</cite>.
 *
 * <p>At the difference of most other map projection providers, this class does not extend {@link MapProjection}
 * because it does not create non-liner kernel. Instead, the projection created by this class is implemented
 * by an affine transform.</p>
 *
 * <p>This provider is <strong>not</strong> suitable for the <cite>Equidistant Cylindrical</cite> projection
 * (EPSG:1028, <span class="deprecated">EPSG:9842</span>). EPSG defines Equidistant Cylindrical projection as
 * the ellipsoidal case of this projection, which uses a more complicated formula than the affine transform
 * used here.</p>
 *
 * <h2>Note</h2>
 * EPSG:1028 and 1029 are the current codes, while EPSG:9842 and 9823 are deprecated codes.
 * The new and deprecated definitions differ only by their parameters. In the Apache SIS implementation,
 * both current and legacy definitions are known, but the legacy names are marked as deprecated.
 *
 * @author  John Grange
 * @author  Martin Desruisseaux (Geomatys)
 *
 * @see PseudoPlateCarree
 * @see <a href="http://geotiff.maptools.org/proj_list/equirectangular.html">GeoTIFF parameters for Equirectangular</a>
 */
@XmlTransient
public final class Equirectangular extends AbstractProvider {
    /**
     * For cross-version compatibility.
     */
    private static final long serialVersionUID = -278288251842178001L;

    /**
     * Name of this projection in EPSG geodetic dataset.
     *
     * @todo Remove with JDK9 after we introduce {@code getInstance()} method.
     */
    public static final String NAME = "Equidistant Cylindrical (Spherical)";

    /*
     * ACCESS POLICY: Only formal EPSG parameters shall be public.
     * Parameters that we add ourselves should be package-privated.
     */

    /**
     * The operation parameter descriptor for the <cite>Latitude of 1st standard parallel</cite> (φ₁) parameter value.
     * Valid values range is (-90 … 90)° and default value is 0°.
     *
     * <!-- Generated by ParameterNameTableGenerator -->
     * <table class="sis">
     *   <caption>Parameter names</caption>
     *   <tr><td> EPSG:    </td><td> Latitude of 1st standard parallel </td></tr>
     *   <tr><td> OGC:     </td><td> standard_parallel_1 </td></tr>
     *   <tr><td> ESRI:    </td><td> Standard_Parallel_1 </td></tr>
     *   <tr><td> NetCDF:  </td><td> standard_parallel </td></tr>
     *   <tr><td> GeoTIFF: </td><td> StdParallel1 </td></tr>
     *   <tr><td> Proj4:   </td><td> lat_ts </td></tr>
     * </table>
     * <b>Notes:</b>
     * <ul>
     *   <li>Value domain: (-90.0 … 90.0)°</li>
     * </ul>
     */
    public static final ParameterDescriptor<Double> STANDARD_PARALLEL;

    /**
     * The operation parameter descriptor for the <cite>Latitude of natural origin</cite> (φ₀) parameter value.
     * Valid values range is (-90 … 90)° and default value is 0°.
     *
     * <p>In theory, this parameter should not be used and its value should be 0 in all cases.
     * This parameter is included for completeness in CRS labeling only, and is declared optional.</p>
     *
     * <!-- Generated by ParameterNameTableGenerator -->
     * <table class="sis">
     *   <caption>Parameter names</caption>
     *   <tr><td> EPSG:    </td><td> Latitude of natural origin </td></tr>
     *   <tr><td> OGC:     </td><td> latitude_of_origin </td></tr>
     *   <tr><td> ESRI:    </td><td> Latitude_Of_Origin </td></tr>
     *   <tr><td> NetCDF:  </td><td> latitude_of_projection_origin </td></tr>
     *   <tr><td> GeoTIFF: </td><td> CenterLat </td></tr>
     *   <tr><td> Proj4:   </td><td> lat_0 </td></tr>
     * </table>
     * <b>Notes:</b>
     * <ul>
     *   <li>Value restricted to 0</li>
     *   <li>Optional</li>
     * </ul>
     */
    static final ParameterDescriptor<Double> LATITUDE_OF_ORIGIN;

    /**
     * The operation parameter descriptor for the <cite>Longitude of natural origin</cite> (λ₀) parameter value.
     * Valid values range is [-180 … 180]° and default value is 0°.
     *
     * <!-- Generated by ParameterNameTableGenerator -->
     * <table class="sis">
     *   <caption>Parameter names</caption>
     *   <tr><td> EPSG:    </td><td> Longitude of natural origin </td></tr>
     *   <tr><td> OGC:     </td><td> central_meridian </td></tr>
     *   <tr><td> ESRI:    </td><td> Central_Meridian </td></tr>
     *   <tr><td> NetCDF:  </td><td> longitude_of_projection_origin </td></tr>
     *   <tr><td> GeoTIFF: </td><td> CenterLong </td></tr>
     *   <tr><td> Proj4:   </td><td> lon_0 </td></tr>
     * </table>
     */
    public static final ParameterDescriptor<Double> LONGITUDE_OF_ORIGIN;

    /**
     * The operation parameter descriptor for the <cite>False easting</cite> (FE) parameter value.
     * Valid values range is unrestricted and default value is 0 metre.
     *
     * <!-- Generated by ParameterNameTableGenerator -->
     * <table class="sis">
     *   <caption>Parameter names</caption>
     *   <tr><td> EPSG:    </td><td> False easting </td></tr>
     *   <tr><td> OGC:     </td><td> false_easting </td></tr>
     *   <tr><td> ESRI:    </td><td> False_Easting </td></tr>
     *   <tr><td> NetCDF:  </td><td> false_easting </td></tr>
     *   <tr><td> GeoTIFF: </td><td> FalseEasting </td></tr>
     *   <tr><td> Proj4:   </td><td> x_0 </td></tr>
     * </table>
     */
    public static final ParameterDescriptor<Double> FALSE_EASTING;

    /**
     * The operation parameter descriptor for the <cite>False northing</cite> (FN) parameter value.
     * Valid values range is unrestricted and default value is 0 metre.
     *
     * <!-- Generated by ParameterNameTableGenerator -->
     * <table class="sis">
     *   <caption>Parameter names</caption>
     *   <tr><td> EPSG:    </td><td> False northing </td></tr>
     *   <tr><td> OGC:     </td><td> false_northing </td></tr>
     *   <tr><td> ESRI:    </td><td> False_Northing </td></tr>
     *   <tr><td> NetCDF:  </td><td> false_northing </td></tr>
     *   <tr><td> GeoTIFF: </td><td> FalseNorthing </td></tr>
     *   <tr><td> Proj4:   </td><td> y_0 </td></tr>
     * </table>
     */
    public static final ParameterDescriptor<Double> FALSE_NORTHING;

    /**
     * The group of all parameters expected by this coordinate operation.
     */
    static final ParameterDescriptorGroup PARAMETERS;
    static {
        final ParameterBuilder builder = builder();

        STANDARD_PARALLEL = createLatitude(builder
                .addIdentifier("8823")
                .addIdentifier(Citations.GEOTIFF, "3078")
                .addName("Latitude of 1st standard parallel")
                .addName(Citations.OGC,     Constants.STANDARD_PARALLEL_1)
                .addName(Citations.ESRI,    "Standard_Parallel_1")
                .addName(Citations.NETCDF,  Constants.STANDARD_PARALLEL)
                .addName(Citations.GEOTIFF, "StdParallel1")
                .addName(Citations.PROJ4,   "lat_ts"), false);

        LONGITUDE_OF_ORIGIN = createLongitude(builder
                .addIdentifier("8802")
                .addIdentifier(Citations.GEOTIFF, "3088")
                .addName("Longitude of natural origin")
                .addName(Citations.OGC,     Constants.CENTRAL_MERIDIAN)
                .addName(Citations.ESRI,    "Central_Meridian")
                .addName(Citations.NETCDF,  "longitude_of_projection_origin")
                .addName(Citations.GEOTIFF, "CenterLong")
                .addName(Citations.PROJ4,   "lon_0"));

        FALSE_EASTING = createShift(builder
                .addIdentifier("8806")
                .addIdentifier(Citations.GEOTIFF, "3082")
                .addName("False easting")
                .addName(Citations.OGC,     Constants.FALSE_EASTING)
                .addName(Citations.ESRI,    "False_Easting")
                .addName(Citations.NETCDF,  Constants.FALSE_EASTING)
                .addName(Citations.GEOTIFF, "FalseEasting")
                .addName(Citations.PROJ4,   "x_0"));

        FALSE_NORTHING = createShift(builder
                .addIdentifier("8807")
                .addIdentifier(Citations.GEOTIFF, "3083")
                .addName("False northing")
                .addName(Citations.OGC,     Constants.FALSE_NORTHING)
                .addName(Citations.ESRI,    "False_Northing")
                .addName(Citations.NETCDF,  Constants.FALSE_NORTHING)
                .addName(Citations.GEOTIFF, "FalseNorthing")
                .addName(Citations.PROJ4,   "y_0"));
        /*
         * "Latitude of natural origin" is not formally parameters of the "Equidistant Cylindrical (Spherical)"
         * projection according EPSG:1029.  But we declare it anyway (as an optional parameter) because it was
         * part of the now deprecated EPSG:9823 definition (and also EPSG:9842, the ellipsoidal case),  and we
         * still see it in use sometimes. However, taking inspiration from the practice done in "Mercator (1SP)"
         * projection, we require that the parameter value must be zero.
         */
        LATITUDE_OF_ORIGIN = createZeroConstant(builder     // Was used by EPSG:9823 (also EPSG:9842).
                .addIdentifier("8801")
                .addIdentifier(Citations.GEOTIFF, "3089")
                .addName("Latitude of natural origin")
                .addName(Citations.OGC,     Constants.LATITUDE_OF_ORIGIN)
                .addName(Citations.ESRI,    "Latitude_Of_Origin")
                .addName(Citations.NETCDF,  "latitude_of_projection_origin")
                .addName(Citations.GEOTIFF, "CenterLat")
                .addName(Citations.PROJ4,   "lat_0")
                .setRemarks(Resources.formatInternational(Resources.Keys.ConstantProjParameterValue_1, 0))
                .setRequired(false));

        // Do not declare the ESRI "Equidistant_Cylindrical" projection name below,
        // for avoiding confusion with EPSG "Equidistant Cylindrical" ellipsoidal projection.
        PARAMETERS = addIdentifierAndLegacy(builder, "1029", "9823")  // 9823 uses deprecated parameter names
                .addName(                   NAME)
                .addName(                   "Plate Carrée")  // Not formally defined by EPSG, but cited in documentation.
                .addName(Citations.OGC,     "Equirectangular")
                .addName(Citations.ESRI,    "Plate_Carree")
                .addName(Citations.GEOTIFF, "CT_Equirectangular")
                .addName(Citations.PROJ4,   "eqc")
                .addIdentifier(Citations.GEOTIFF, "17")
                .createGroupForMapProjection(
                        STANDARD_PARALLEL,
                        LATITUDE_OF_ORIGIN,     // Not formally an Equirectangular parameter.
                        LONGITUDE_OF_ORIGIN,
                        FALSE_EASTING,
                        FALSE_NORTHING);
    }

    /**
     * Constructs a new provider.
     *
     * @see MapProjection#MapProjection(Class, ParameterDescriptorGroup)
     */
    public Equirectangular() {
        super(CylindricalProjection.class, PARAMETERS,
              EllipsoidalCS.class, true,
              CartesianCS.class,   false);
    }

    /**
     * Gets a parameter value identified by the given descriptor and stores it only if different than zero.
     * This method performs the same work than {@code Initializer.getAndStore(ParameterDescriptor)} in the
     * {@link org.apache.sis.referencing.operation.projection} package.
     *
     * @param  source      the parameters from which to read the value.
     * @param  target      where to store the parameter values.
     * @param  descriptor  the descriptor that specify the parameter names and desired units.
     * @return the parameter value in the units given by the descriptor.
     * @throws IllegalArgumentException if the given value is out of bounds.
     */
    private static double getAndStore(final Parameters source, final ParameterValueGroup target,
            final ParameterDescriptor<Double> descriptor) throws IllegalArgumentException
    {
        final double value = source.doubleValue(descriptor);    // Apply a unit conversion if needed.
        MapProjection.validate(descriptor, value);              // Unconditional validation for all parameters.
        if (value != 0) {                                       // All default values in this class are zero.
            target.parameter(descriptor.getName().getCode()).setValue(value);
        }
        return value;
    }

    /**
     * Creates an Equirectangular projection from the specified group of parameter values. This method is an
     * adaptation of {@link org.apache.sis.referencing.operation.projection.NormalizedProjection} constructor,
     * reproduced in this method because we will create an affine transform instead of the usual projection
     * classes.
     *
     * @param  factory     the factory to use if this constructor needs to create other math transforms.
     * @param  parameters  the parameter values that define the transform to create.
     * @return the map projection created from the given parameter values.
     * @throws FactoryException if an error occurred while creating the math transform.
     */
    @Override
    public MathTransform createMathTransform(final MathTransformFactory factory, final ParameterValueGroup parameters)
            throws FactoryException
    {
        final Parameters p = Parameters.castOrWrap(parameters);
        final ContextualParameters context = new ContextualParameters(PARAMETERS, 2, 2);
        double a  = getAndStore(p, context, MapProjection.SEMI_MAJOR);
        double b  = getAndStore(p, context, MapProjection.SEMI_MINOR);
        double λ0 = getAndStore(p, context, LONGITUDE_OF_ORIGIN);
        double φ0 = getAndStore(p, context, LATITUDE_OF_ORIGIN);
        double φ1 = getAndStore(p, context, STANDARD_PARALLEL);
        double fe = getAndStore(p, context, FALSE_EASTING);
        double fn = getAndStore(p, context, FALSE_NORTHING);
        /*
         * Perform following transformation, in that order. Note that following
         * AffineTransform convention, the Java code appears in reverse order:
         *
         *   1) Subtract φ0 to the latitude.
         *   2) Subtract λ0 to the longitude.
         *   3) Convert degrees to radians.
         *   4) Scale longitude by cos(φ1).
         */
        φ1 = toRadians(φ1);
        final MatrixSIS normalize = context.getMatrix(ContextualParameters.MatrixRole.NORMALIZATION);
        normalize.convertBefore(0, cos(φ1), null);
        context.normalizeGeographicInputs(λ0)
               .convertBefore(1, null, -φ0);
        /*
         * At this point, we usually invoke 'denormalize.convertAfter(…, a, …)' where 'a' (the semi-major axis length)
         * is taken as the Earth radius (R). However, quoting EPSG: "If the figure of the earth used is an ellipsoid
         * rather than a sphere then R should be calculated as the radius of the conformal sphere at the projection
         * origin at latitude φ1 using the formula for RC given in section 1.2, table 3".
         */
        if (a != b) {
            final double rs = b / a;
            final double sinφ1 = sin(φ1);
            a = b / (1 - (1 - rs*rs) * (sinφ1*sinφ1));
        }
        final DoubleDouble k = DoubleDouble.of(a, true);
        final MatrixSIS denormalize = context.getMatrix(ContextualParameters.MatrixRole.DENORMALIZATION);
        denormalize.convertAfter(0, k, DoubleDouble.of(fe, true));
        denormalize.convertAfter(1, k, DoubleDouble.of(fn, true));
        /*
         * Creates the ConcatenatedTransform, letting the factory returns the cached instance
         * if the caller already invoked this method previously (which usually do not happen).
         */
        MathTransform mt = context.completeTransform(factory, MathTransforms.identity(2));
        if (mt instanceof AffineTransform) {  // Always true in Apache SIS implementation.
            mt = new ParameterizedAffine((AffineTransform) mt, context, true);
        }
        return mt;
    }
}
