using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace tileserver
{
    class Mercator
    {
        /**
        * Tile size.<br>
        */
        public static int TILE_SIZE = 256;

        /**
        * Maximum latitude.<br>
        */
        private static double MAX_LAT = 85.05112877980659;

        /**
        * Minimum latitude.<br>
        */
        private static double MIN_LAT = -85.05112877980659;
    
        /**
        * Equatorial earth radius for EPSG:3857 (Mercator).<br>
        */
        private static double EARTH_RADIUS = 6378137;


        /**
        * Returns earth radius size at zoom level.<br>
        * @param aZoomlevel Zoom level.
        * @return Earth radius size at zoom level.
        */

        public static double radius(int aZoomlevel)
        {
            return (TILE_SIZE * (1 << aZoomlevel)) / (2.0 * Math.PI);
        }




        /**
        * Returns the absolut number of pixels in y or x, defined as: 2^Zoomlevel * TILE_WIDTH where TILE_WIDTH is the width of a tile in pixels.<br>
        * @param aZoomlevel zoom level to request pixel data.
        * @return number of pixels.
        */

        public static int getMaxPixels(int aZoomlevel)
        {
            return TILE_SIZE * (1 << aZoomlevel);
        }


        public static int falseEasting(int aZoomlevel)
        {
            return getMaxPixels(aZoomlevel) / 2;
        }

        public static int falseNorthing(int aZoomlevel)
        {
            return (-1 * getMaxPixels(aZoomlevel) / 2);
        }



        /**
        * Transform pixelspace to coordinates and get the distance.
        *
        * @param x1 the first x coordinate.
        * @param y1 the first y coordinate.
        * @param x2 the second x coordinate.
        * @param y2 the second y coordinate.
        * @param zoomLevel the zoom level.
        * @return the distance.
        */


        public static double getDistance(int x1, int y1, int x2, int y2, int zoomLevel)
        {
            double la1 = YToLatitude(y1, zoomLevel);
            double lo1 = XToLongitude(x1, zoomLevel);
            double la2 = YToLatitude(y2, zoomLevel);
            double lo2 = XToLongitude(x2, zoomLevel);

            return getDistance(la1, lo1, la2, lo2);
        }





        /**
        * Gets the distance using Spherical law of cosines.
        *
        * @param la1 the Latitude in degrees.
        * @param lo1 the Longitude in degrees.
        * @param la2 the Latitude from 2nd coordinate in degrees.
        * @param lo2 the Longitude from 2nd coordinate in degrees.
        * @return the distance.
        * 
        */

        public static double getDistance(double la1, double lo1, double la2, double lo2)
        {
            double aStartLat = ToRadians(la1);
            double aStartLong = ToRadians(lo1);
            double aEndLat = ToRadians(la2);
            double aEndLong = ToRadians(lo2);
            double distance = Math.Acos(Math.Sin(aStartLat) * Math.Sin(aEndLat) + Math.Cos(aStartLat) * Math.Cos(aEndLat) * Math.Cos(aEndLong - aStartLong));
            return (EARTH_RADIUS * distance);
        }



        /**
        * Transform longitude to pixelspace
        *
        * <p>
        * Mathematical optimization<br>
        * <code>
        * x = radius(aZoomlevel) * toRadians(aLongitude) + falseEasting(aZoomLevel)<br>
        * x = getMaxPixels(aZoomlevel) / (2 * PI) * (aLongitude * PI) / 180 + getMaxPixels(aZoomlevel) / 2<br>
        * x = getMaxPixels(aZoomlevel) * aLongitude / 360 + 180 * getMaxPixels(aZoomlevel) / 360<br>
        * x = getMaxPixels(aZoomlevel) * (aLongitude + 180) / 360<br>
        * </code>
        * </p>
        *
        * @param aLongitude [-180..180]
        * @return [0..2^Zoomlevel*TILE_SIZE]
        */
    
        public static double longitudeToX(double aLongitude, int aZoomlevel) 
        {
            int mp = getMaxPixels(aZoomlevel);
            double x = (mp * (aLongitude + 180L)) / 360L;
            return Math.Min(x, mp - 1);
        }

    
    

    
        /**
        * Transform longitude to pixelspace
        *
        * <p>
        * Mathematical optimization<br>
        * <code>
        * x = radius(aZoomlevel) * toRadians(aLongitude) + falseEasting(aZoomLevel)<br>
        * x = getMaxPixels(aZoomlevel) / (2 * PI) * (aLongitude * PI) / 180 + getMaxPixels(aZoomlevel) / 2<br>
        * x = getMaxPixels(aZoomlevel) * aLongitude / 360 + 180 * getMaxPixels(aZoomlevel) / 360<br>
        * x = getMaxPixels(aZoomlevel) * (aLongitude + 180) / 360<br>
        * </code>
        * </p>
        *
        * @param aLongitude [-180..180]
        * @return [0..2^Zoomlevel]
        */
    
    
        public static int longitudeToTileX(double aLongitude, int aZoomlevel) 
        {
            return (int)(longitudeToX(aLongitude, aZoomlevel) / (double)TILE_SIZE);
        }
    
    
        /**
        * Transforms latitude to pixelspace
        * <p>
        * Mathematical optimization<br>
        * <code>
        * log(u) := log((1.0 + sin(toRadians(aLat))) / (1.0 - sin(toRadians(aLat))<br>
        *
        * y = -1 * (radius(aZoomlevel) / 2 * log(u)))) - falseNorthing(aZoomlevel))<br>
        * y = -1 * (getMaxPixel(aZoomlevel) / 2 * PI / 2 * log(u)) - -1 * getMaxPixel(aZoomLevel) / 2<br>
        * y = getMaxPixel(aZoomlevel) / (-4 * PI) * log(u)) + getMaxPixel(aZoomLevel) / 2<br>
        * y = getMaxPixel(aZoomlevel) * ((log(u) / (-4 * PI)) + 1/2)<br>
        * </code>
        * </p>
        * @param aLat [-90...90]
        * @return [0..2^Zoomlevel*TILE_SIZE]
        */
    
        public static double latitudeToY(double aLat, int aZoomlevel) 
        {
            if (aLat < MIN_LAT)
            {
                aLat = MIN_LAT;
            }
            else if (aLat > MAX_LAT)
            {
                aLat = MAX_LAT;
            }
            double sinLat = Math.Sin(ToRadians(aLat));
            double log = Math.Log((1.0 + sinLat) / (1.0 - sinLat));
            int mp = getMaxPixels(aZoomlevel);
            double y = mp * (0.5 - (log / (4.0 * Math.PI)));
            return Math.Min(y, mp - 1);
        }

    
        /**
        * Transforms latitude to pixelspace
        * <p>
        * Mathematical optimization<br>
        * <code>
        * log(u) := log((1.0 + sin(toRadians(aLat))) / (1.0 - sin(toRadians(aLat))<br>
        *
        * y = -1 * (radius(aZoomlevel) / 2 * log(u)))) - falseNorthing(aZoomlevel))<br>
        * y = -1 * (getMaxPixel(aZoomlevel) / 2 * PI / 2 * log(u)) - -1 * getMaxPixel(aZoomLevel) / 2<br>
        * y = getMaxPixel(aZoomlevel) / (-4 * PI) * log(u)) + getMaxPixel(aZoomLevel) / 2<br>
        * y = getMaxPixel(aZoomlevel) * ((log(u) / (-4 * PI)) + 1/2)<br>
        * </code>
        * </p>
        * @param aLat [-90...90]
        * @return [0..2^Zoomlevel]
        */
    
        public static int latitudeToTileY(double aLat, int aZoomlevel) 
        {
            return (int)(latitudeToY(aLat, aZoomlevel) / TILE_SIZE);
        }
    
    
    
        /**
        * Transforms pixel coordinate X to longitude
        *
        * <p>
        * Mathematical optimization<br>
        * <code>
        * lon = toDegree((aX - falseEasting(aZoomlevel)) / radius(aZoomlevel))<br>
        * lon = 180 / PI * ((aX - getMaxPixels(aZoomlevel) / 2) / getMaxPixels(aZoomlevel) / (2 * PI)<br>
        * lon = 180 * ((aX - getMaxPixels(aZoomlevel) / 2) / getMaxPixels(aZoomlevel))<br>
        * lon = 360 / getMaxPixels(aZoomlevel) * (aX - getMaxPixels(aZoomlevel) / 2)<br>
        * lon = 360 * aX / getMaxPixels(aZoomlevel) - 180<br>
        * </code>
        * </p>
        * @param aX [0..2^Zoomlevel*TILE_WIDTH]
        * @return [-180..180]
        */
    
        public static double XToLongitude(int aX, int aZoomlevel) 
        {
            return ((360d * aX) / getMaxPixels(aZoomlevel)) - 180.0;
        }

    
    
        /**
        * Transforms pixel coordinate Y to latitude
        *
        * @param aY [0..2^Zoomlevel*TILE_WIDTH]
        * @return [MIN_LAT..MAX_LAT] is about [-85..85]
        */
    
        public static double YToLatitude(int aY, int aZoomlevel) 
        {
            aY += falseNorthing(aZoomlevel);
            double latitude = (Math.PI / 2) - (2 * Math.Atan(Math.Exp(-1.0 * aY / radius(aZoomlevel))));
            return -1 * ToDegrees(latitude);
        }

    
 
    
        /**
         * Ecuaciones para los metros por cada píxel
         * La distancia representada por un píxel (S) está dada por
         * 
         * S=C*cos(y)/2^(z+8)
         * *
         * donde...
         * C es la circunferencia ecuatorial de la Tierra
         * z es el nivel de zoom
         * y es la latitud del lugar cuya escala cartográfica queremos saber.
         * 
         * Asegúrese de que su calculadora esté en modo sexagesimal para los ángulos, 
         * a menos que por algún motivo desee expresar la latitud en radianes. 
         * C debe estar expresado en cualquier unidad de su interés (millas, metros, pies, cualquiera). 
         * Puesto que la Tierra tiene en realidad forma elipsoidal, habrá un pequeño error en este cálculo. 
         * Pero es muy pequeño. (0.3% es el error máximo) 
         * 
         *	grado	m / pixel	escala	
         *  1				
         *  2				
         *  3				
         *  4				
         *  5				
         *  6				
         *  7				
         *  8	1,406	610,984	1 :	2.000.000
         *  9	0,703	305,492	1 :	1.100.000
         *  10	0,352	152,746	1 :	500.000
         *  11	0,176	76,373	1 :	250.000
         *  12	0,088	38,187	1 :	150.000
         *  13	0,044	19,093	1 :	70.000
         *  14	0,022	9,547	1 :	35.000
         *  15	0,011	4,773	1 :	15.000
         *  16	0,005	2,387	1 :	8.000
         *  17	0,003	1,193	1 :	4.000
         *  18	0,001	0,597	1 :	2.000
         * 
         */
    
        public static double distanceOnePixel(double latitude, int zoom)
        {
            double C = Math.PI * 2 * EARTH_RADIUS;
            double S = (C * Math.Cos(latitude) / Math.Pow( 2, zoom + 8));
            return S;
        }



        public static double distanceOnePixel(int x, int y, int zoom)
        {
            double lon1 = XToLongitude(x, zoom);
            double lat1 = YToLatitude(y, zoom);
            double lon2 = XToLongitude(x + 1, zoom);
            double lat2 = YToLatitude(y, zoom);
            double distance = getDistance(lat1, lon1, lat2, lon2);
            return distance;
        }



        private static double ToRadians(double degrees)
        {
            return (Math.PI / 180) * degrees;
        }



        private static double ToDegrees(double radians)
        {
            return (180.0 / Math.PI) * radians;
        }

    

    }
}