<?php
declare(strict_types=1);

namespace algorithm;

/**
 * 地理算法
 */
class Geometry
{

    
    /**
     * 火星坐标Gcj02 转 WGS84
     * 
     * @param double $lat 纬度
     * @param double $lng 经度
     */
    public static function gcj02ToWgs84($lat, $lng)
    {
       
        $pi = 3.14159265358979324;
        $a = 6378245.0;
        $ee = 0.00669342162296594323;
        
        $dlat = self::wgsTOgcjTransformLat($lng - 105.0, $lat - 35.0);
        $dlng = self::wgsTOgcjTransformLon($lng - 105.0, $lat - 35.0);
        $radlat = $lat / 180.0 * $pi;
        $magic = sin($radlat);
        $magic = 1 - $ee * $magic * $magic;
        $sqrtmagic = sqrt($magic);
        $dlat = ($dlat * 180.0) / (($a * (1 - $ee)) / ($magic * $sqrtmagic) * $pi);
        $dlng = ($dlng * 180.0) / ($a / $sqrtmagic * cos($radlat) * $pi);
        $mglat = $lat + $dlat;
        $mglng = $lng + $dlng;
        return [ 
            'lng' =>round($lng * 2 - $mglng,6),
            'lat' =>round($lat * 2 - $mglat,6), 
        ];
    }

    /**
     * 火星坐标Gcj02 转 百度地图坐标BD09
     * 腾讯地图用的也是GCJ02坐标
     * @param double $lat 纬度
     * @param double $lng 经度
     */
    public static function convertGcj02ToBd09($lat, $lng)
    {
        $lat = floatval($lat);
        $lng = floatval($lng);
        $x_pi = 3.14159265358979324 * 3000.0 / 180.0;
        $x = $lng;
        $y = $lat;
        $z =sqrt($x * $x + $y * $y) + 0.00002 * sin($y * $x_pi);
        $theta = atan2($y, $x) + 0.000003 * cos($x * $x_pi);
        $lng = $z * cos($theta) + 0.0065;
        $lat = $z * sin($theta) + 0.006;
        $lng = round($lng,6);
        $lat = round($lat,6);
        return array('lng'=>$lng,'lat'=>$lat);
    }

    /**
     * 百度地图坐标Bd09 转 火星坐标Gcj02
     * @param double $lat 纬度
     * @param double $lng 经度
     */
    public static function convertBd09ToGcj02($lat, $lng)
    {
        $lat = floatval($lat);
        $lng = floatval($lng);
        $x_pi = 3.14159265358979324 * 3000.0 / 180.0;
        $x = $lng - 0.0065;
        $y = $lat - 0.006;
        $z = sqrt($x * $x + $y * $y) - 0.00002 * sin($y * $x_pi);
        $theta = atan2($y, $x) - 0.000003 * cos($x * $x_pi);
        $lng = $z * cos($theta);
        $lat = $z * sin($theta);
        $lng = round($lng,6);
        $lat = round($lat,6);
        return array('lng'=>$lng,'lat'=>$lat);
    }

    /**
     * wgs84坐标转国测局坐标
     * @param double $lat 纬度
     * @param double $lng 经度
     */
    public static function wgsTOgcj($lat, $lng)
    {
        $pi = 3.14159265358979324;
        $a = 6378245.0;
        $ee = 0.00669342162296594323;
        $wgLat = $lat;
        $wgLon = $lng;
        if (self::outOfChina($wgLat, $wgLon)){
            return array('lat'=>$wgLat, 'lng'=>$wgLon);
        }
        $dLat = self::wgsTOgcjTransformLat($wgLon - 105.0, $wgLat - 35.0);
        $dLon = self::wgsTOgcjTransformLon($wgLon - 105.0, $wgLat - 35.0);
        $radLat = $wgLat / 180.0 * $pi;
        $magic = sin($radLat);
        $magic = 1 - $ee * $magic * $magic;
        $sqrtMagic = sqrt($magic);
        $dLat = ($dLat * 180.0) / (($a * (1 - $ee)) / ($magic * $sqrtMagic) * $pi);
        $dLon = ($dLon * 180.0) / ($a / $sqrtMagic * cos($radLat) * $pi);

        $mgLat = $wgLat + $dLat;
        $mgLon = $wgLon + $dLon;

        return array('lat'=>$mgLat, 'lng'=>$mgLon);
    }

    /**
     * 经纬度是否超出中国范围
     * @param double $lat 纬度
     * @param double $lng 经度
     */
    public static function outOfChina($lat, $lng)
    {
        if ($lng < 72.004 || $lng > 137.8347)
            return true;
        if ($lat < 0.8293 || $lat > 55.8271)
            return true;

        return false;
    }
    
     
  
    
    /**
     * wgsTOgcj转换纬度
     * @param double $lat 纬度
     * @param double $lng 经度
     */
    private static function wgsTOgcjTransformLat( $x, $y ){
        $pi = 3.14159265358979324;
        $ret = -100.0 + 2.0 * $x + 3.0 * $y + 0.2 * $y * $y + 0.1 * $x * $y + 0.2 * sqrt(abs($x));
        $ret += (20.0 * sin(6.0 * $x * $pi) + 20.0 * sin(2.0 * $x * $pi)) * 2.0 / 3.0;
        $ret += (20.0 * sin($y * $pi) + 40.0 * sin($y / 3.0 * $pi)) * 2.0 / 3.0;
        $ret += (160.0 * sin($y / 12.0 * $pi) + 320 * sin($y * $pi / 30.0)) * 2.0 / 3.0;
        return $ret;
    }

    /**
     * wgsTOgcj转换经度
     * @param double $lat 纬度
     * @param double $lng 经度
     */
    private static function wgsTOgcjTransformLon( $x, $y ){
        $pi = 3.14159265358979324;
        $ret = 300.0 + $x + 2.0 * $y + 0.1 * $x * $x + 0.1 * $x * $y + 0.1 * sqrt(abs($x));
        $ret += (20.0 * sin(6.0 * $x * $pi) + 20.0 * sin(2.0 * $x * $pi)) * 2.0 / 3.0;
        $ret += (20.0 * sin($x * $pi) + 40.0 * sin($x / 3.0 * $pi)) * 2.0 / 3.0;
        $ret += (150.0 * sin($x / 12.0 * $pi) + 300.0 * sin($x / 30.0 * $pi)) * 2.0 / 3.0;
        return $ret;
    }

     /**
     * 获取百度2点间距离(单位：km)
     * @param  $point1 坐标1
     * @param  $point2 坐标2
     */
    public static function distanceBetween2BdPoints($point1, $point2){
        $PI = 3.1415926535898;
        $earthRadius = 6371.393;//单位：km
        $radLat1 = $point1['lat'] * ($PI / 180);
        $radLat2 = $point2['lat'] * ($PI / 180);
        $a = $radLat1 - $radLat2;
        $b = ($point1['lng'] * ($PI / 180)) - ($point2['lng'] * ($PI / 180));
        $s = 2 * asin(sqrt(pow(sin($a/2),2) + cos($radLat1)*cos($radLat2)*pow(sin($b/2),2)));
        $s = $s * $earthRadius;
        $s = round($s * 10000) / 10000;
        return $s;
    }

    /**
     *坐标点是否在多边形范围内
     * @param  $polygon 多边形点数组
     * @param  $lnglat 经纬度数组
     */
    public static function isInPolygon($polygon, $lnglat){
        $count = count($polygon);
        $px = $lnglat['lat'];
        $py = $lnglat['lng'];
        $flag = FALSE;
        for ($i = 0, $j = $count - 1; $i < $count; $j = $i, $i++) {
            $sy = $polygon[$i]['lng'];
            $sx = $polygon[$i]['lat'];
            $ty = $polygon[$j]['lng'];
            $tx = $polygon[$j]['lat'];
            if ($px == $sx && $py == $sy || $px == $tx && $py == $ty)
                return TRUE;
            if ($sy < $py && $ty >= $py || $sy >= $py && $ty < $py) {
                $x = $sx + ($py - $sy) * ($tx - $sx) / ($ty - $sy); 
                if ($x == $px) return TRUE; 
                if ($x > $px)
                    $flag = !$flag;
            }
        }
        return $flag;
    }
}