Chinese GPS coordinate calculate

GPS.kt

import java.math.BigDecimal
import kotlin.math.*

/**
 * No need to apply transformations on non-Chinese coordinates.
 *
 * This is not follow the geo definition line. It just thinks a square.
 * */
fun isOutOfChina(lat: Double, lng: Double): Boolean {
    if (lng < 72.004 || lng > 137.8347) {
        return true
    }
    return lat < 0.8293 || lat > 55.8271
}

// credit: https://github.com/googollee/eviltransform
private fun transform(x: Double, y: Double): Pair<Double, Double> {
    val xy = x * y
    val absX = sqrt(abs(x))
    val xPi = x * PI
    val yPi = y * PI
    val d = 20.0 * sin(6.0 * xPi) + 20.0 * sin(2.0 * xPi)

    var resultLat = d
    var resultLng = d

    resultLat += 20.0 * sin(yPi) + 40.0 * sin(yPi / 3.0)
    resultLng += 20.0 * sin(xPi) + 40.0 * sin(xPi / 3.0)

    resultLat += 160.0 * sin(yPi / 12.0) + 320 * sin(yPi / 30.0)
    resultLng += 150.0 * sin(xPi / 12.0) + 300.0 * sin(xPi / 30.0)

    resultLat *= 2.0 / 3.0
    resultLng *= 2.0 / 3.0

    resultLat += -100.0 + 2.0 * x + 3.0 * y + 0.2 * y * y + 0.1 * xy + 0.2 * absX
    resultLng += 300.0 + x + 2.0 * y + 0.1 * x * x + 0.1 * xy + 0.1 * absX

    return resultLat to resultLng
}

// credit: https://github.com/googollee/eviltransform
private fun delta(lat: Double, lng: Double): Pair<Double, Double> {
    val ee = BigDecimal("0.00669342162296594323")
    var (dLat, dLng) = transform(lng - 105.0, lat - 35.0)
    val radLat = lat / 180.0 * PI
    var magic = sin(radLat)
    magic = 1.0 - (ee * magic.toBigDecimal() * magic.toBigDecimal()).toDouble()
    val sqrtMagic = sqrt(magic)
    // someone suggests the earthR should be this, but we use the original value here
//    val earthR = 6378245.0.toBigDecimal()
    val earthR = 6378137.0.toBigDecimal()
    dLat = (dLat * 180.0) / ((earthR * (BigDecimal.ONE - ee)).toDouble() / (magic * sqrtMagic) * PI)
    dLng = (dLng * 180.0) / (earthR.toDouble() / sqrtMagic * cos(radLat) * PI)

    return dLat to dLng
}

/**
 * Convert wgs coordinate (standard) to gcj coordinate (China).
 * */
fun wgs2gcj(wgsLat: Double, wgsLng: Double): Pair<Double, Double> {
    if (isOutOfChina(wgsLat, wgsLng)) return wgsLat to wgsLng
    val (dLat, dLng) = delta(wgsLat, wgsLng)
    return (wgsLat + dLat) to (wgsLng + dLng)
}

/**
 * Convert gcj coordinate (China) to wgs coordinate (standard).
 * */
fun gcj2wgs(gcjLat: Double, gcjLng: Double): Pair<Double, Double> {
    if (isOutOfChina(gcjLat, gcjLng)) return gcjLat to gcjLng
    val (dLat, dLng) = delta(gcjLat, gcjLng)
    return (gcjLat - dLat) to (gcjLng - dLng)
}

/**
 * Here is a more precise implementation.
 * */
fun gcj2wgsPrecise(gcjLat: Double, gcjLng: Double): Pair<Double, Double> {
    val initialDelta = 0.01
    val threshold = 0.000001
    var dLat = initialDelta
    var dLng = initialDelta
    var mLat = gcjLat - dLat
    var mLng = gcjLng - dLng
    var pLat = gcjLat + dLat
    var pLng = gcjLng + dLng

    var wgsLat = 0.0
    var wgsLng = 0.0

    for (i in 0 until 30) {
        wgsLat = (mLat + pLat) / 2
        wgsLng = (mLng + pLng) / 2

        val (tmpLat, tmpLng) = wgs2gcj(wgsLat, wgsLng)
        dLat = tmpLat - gcjLat
        dLng = tmpLng - gcjLng

        if (abs(dLat) < threshold && abs(dLng) < threshold)
            break

        if (dLat > 0) pLat = wgsLat else mLat = wgsLat
        if (dLng > 0) pLng = wgsLng else mLng = wgsLng
    }

    return wgsLat to wgsLng
}

/**
 * Convert GCJ coordinate (China) to bd coordinate (Baidu).
 * */
fun gcj2bd(gcjLat: Double, gcjLng: Double): Pair<Double, Double> {
    val xPi = PI * 3000.0 / 180.0
    val z = sqrt(gcjLng * gcjLng + gcjLat * gcjLat) + 0.00002 * sin(gcjLat * xPi)
    val theta = atan2(gcjLat, gcjLng) + 0.000003 * cos(gcjLng * xPi)
    return (z * sin(theta) + 0.006) to (z * cos(theta) + 0.0065)
}

// test
fun main() {
    val lat = 39.909183
    val lng = 116.397458
    println(gcj2wgs(lat, lng))
    println(gcj2wgsPrecise(lat, lng))
    println(wgs2gcj(lat, lng))
    println(wgs2gcj(lat, lng).let { gcj2bd(it.first, it.second) })
}