diff --git a/docs/circle/CircleCompose.kt.md b/docs/circle/CircleCompose.kt.md index 2eacaa2..2ed96dc 100644 --- a/docs/circle/CircleCompose.kt.md +++ b/docs/circle/CircleCompose.kt.md @@ -14,7 +14,7 @@ fun MapViewScope.Circle( radiusMeters: Double, id: String? = null, strokeColor: Color = Color.Red, - strokeWidth: Dp = 2.dp, + strokeWidth: Dp = 1.dp, fillColor: Color = Color.White.copy(alpha = 0.5f), zIndex: Int? = null, extra: Serializable? = null, diff --git a/docs/groundimage/GroundImageComponent.kt.md b/docs/groundimage/GroundImageComponent.kt.md index 4b00049..f24a2a6 100644 --- a/docs/groundimage/GroundImageComponent.kt.md +++ b/docs/groundimage/GroundImageComponent.kt.md @@ -21,7 +21,7 @@ its geographical bounds, and other optional properties directly as parameters. fun MapViewScope.GroundImage( bounds: GeoRectBounds, image: Drawable, - opacity: Float = 0.5f, + opacity: Float = 1.0f, tileSize: Int = GroundImageTileProvider.DEFAULT_TILE_SIZE, id: String? = null, extra: Serializable? = null, @@ -46,7 +46,7 @@ be called within the `MapView` composable's content scope. - `opacity` - Type: `Float` - Description: The opacity of the image, ranging from `0.0` (fully transparent) to `1.0` (fully - opaque). Defaults to `0.5f`. + opaque). Defaults to `1.0f`. - `tileSize` - Type: `Int` - Description: The size of the tiles (in pixels) used to render the image on the map. Defaults diff --git a/docs/polygon/Polygon.kt.md b/docs/polygon/Polygon.kt.md index df2628f..06670cc 100644 --- a/docs/polygon/Polygon.kt.md +++ b/docs/polygon/Polygon.kt.md @@ -42,7 +42,7 @@ The `PolygonState` is initialized with the following parameters: - Description: The color of the polygon's outline. - `strokeWidth` - Type: `Dp` - - Default: `2.dp` + - Default: `1.dp` - Description: The width of the polygon's outline. - `fillColor` - Type: `Color` diff --git a/docs/spherical/Spherical.kt.md b/docs/spherical/Spherical.kt.md index 2318e8c..8df2686 100644 --- a/docs/spherical/Spherical.kt.md +++ b/docs/spherical/Spherical.kt.md @@ -313,8 +313,11 @@ fun computeSignedArea(path: List): Double This function computes the signed area of a closed polygon. The sign of the result indicates the orientation of the vertices: -- **Positive value:** The vertices are in a counter-clockwise order. -- **Negative value:** The vertices are in a clockwise order. +- **Positive value:** The vertices are in a clockwise order. +- **Negative value:** The vertices are in a counter-clockwise order. + +Orientation is as seen on a north-up map (longitude increasing to the right, latitude increasing +upward). The magnitude of the result is the area in square meters. The path is assumed to be closed. @@ -327,28 +330,28 @@ The magnitude of the result is the area in square meters. The path is assumed to ### Returns - Type: `Double` -- Description: The signed area in square meters. Positive for CCW, negative for CW. +- Description: The signed area in square meters. Positive for CW, negative for CCW. ### Example ```kotlin -// Counter-clockwise polygon -val ccwPolygon = listOf( +// Clockwise polygon (SW -> NW -> NE -> SE) +val cwPolygon = listOf( GeoPoint(latitude = 40.0, longitude = -75.0), GeoPoint(latitude = 41.0, longitude = -75.0), GeoPoint(latitude = 41.0, longitude = -74.0), GeoPoint(latitude = 40.0, longitude = -74.0) ) -val signedAreaCCW = Spherical.computeSignedArea(ccwPolygon) -// signedAreaCCW will be a positive value -println("Counter-clockwise signed area: $signedAreaCCW") - -// Clockwise polygon -val cwPolygon = ccwPolygon.reversed() val signedAreaCW = Spherical.computeSignedArea(cwPolygon) -// signedAreaCW will be a negative value with the same magnitude +// signedAreaCW will be a positive value println("Clockwise signed area: $signedAreaCW") + +// Counter-clockwise polygon +val ccwPolygon = cwPolygon.reversed() +val signedAreaCCW = Spherical.computeSignedArea(ccwPolygon) +// signedAreaCCW will be a negative value with the same magnitude +println("Counter-clockwise signed area: $signedAreaCCW") ``` --- diff --git a/src/main/java/com/mapconductor/core/polygon/Polygon.kt b/src/main/java/com/mapconductor/core/polygon/Polygon.kt index ea15879..685bbba 100644 --- a/src/main/java/com/mapconductor/core/polygon/Polygon.kt +++ b/src/main/java/com/mapconductor/core/polygon/Polygon.kt @@ -15,7 +15,7 @@ class PolygonState( holes: List> = emptyList(), id: String? = null, strokeColor: Color = Color.Black, - strokeWidth: Dp = 2.dp, + strokeWidth: Dp = 1.dp, fillColor: Color = Color.Transparent, geodesic: Boolean = false, clickable: Boolean = true, diff --git a/src/main/java/com/mapconductor/core/spherical/Spherical.kt b/src/main/java/com/mapconductor/core/spherical/Spherical.kt index b27bcb8..efeac65 100644 --- a/src/main/java/com/mapconductor/core/spherical/Spherical.kt +++ b/src/main/java/com/mapconductor/core/spherical/Spherical.kt @@ -1,286 +1,287 @@ -package com.mapconductor.core.spherical - -import com.mapconductor.core.features.GeoPoint -import com.mapconductor.core.features.GeoPointInterface -import com.mapconductor.core.projection.Earth -import kotlin.math.abs -import kotlin.math.acos -import kotlin.math.asin -import kotlin.math.atan2 -import kotlin.math.cos -import kotlin.math.sin -import kotlin.math.sqrt - -/** - * Spherical geometry utility functions for calculating distances, headings, - * and positions on Earth's surface using the spherical Earth model. - * - * This is a Kotlin port of the Cordova Google Maps plugin spherical.js utilities. - * Uses GeoPointInterface instead of LatLng for coordinate representation. - */ -object Spherical { - // Mathematical constants - private const val PI = Math.PI - private const val RAD_TO_DEG = 180.0 / PI - private const val DEG_TO_RAD = PI / 180.0 - - /** - * Returns the distance, in meters, between two GeoPointInterface locations. - * Uses the haversine formula. - * - * @param from Starting point - * @param to Ending point - * @return Distance in meters - */ - fun computeDistanceBetween( - from: GeoPointInterface, - to: GeoPointInterface, - ): Double { - val lat1Rad = from.latitude * DEG_TO_RAD - val lat2Rad = to.latitude * DEG_TO_RAD - val deltaLat = (to.latitude - from.latitude) * DEG_TO_RAD - val deltaLng = (to.longitude - from.longitude) * DEG_TO_RAD - - val haversineA = - sin(deltaLat / 2) * sin(deltaLat / 2) + - cos(lat1Rad) * cos(lat2Rad) * - sin(deltaLng / 2) * sin(deltaLng / 2) - - val centralAngle = 2 * atan2(sqrt(haversineA), sqrt(1 - haversineA)) - - return Earth.RADIUS_METERS * centralAngle - } - - /** - * Returns the heading from one GeoPointInterface to another GeoPointInterface. - * Headings are expressed in degrees clockwise from North within the range (-180, 180). - * - * @param from Starting point - * @param to Ending point - * @return Heading in degrees - */ - fun computeHeading( - from: GeoPointInterface, - to: GeoPointInterface, - ): Double { - val lat1Rad = from.latitude * DEG_TO_RAD - val lat2Rad = to.latitude * DEG_TO_RAD - val deltaLng = (to.longitude - from.longitude) * DEG_TO_RAD - - val deltaY = sin(deltaLng) * cos(lat2Rad) - val deltaX = cos(lat1Rad) * sin(lat2Rad) - sin(lat1Rad) * cos(lat2Rad) * cos(deltaLng) - - var heading = atan2(deltaY, deltaX) * RAD_TO_DEG - - // Normalize to (-180, 180] - while (heading > 180) heading -= 360 - while (heading <= -180) heading += 360 - - return heading - } - - /** - * Returns the GeoPointInterface resulting from moving a distance from an origin - * in the specified heading (expressed in degrees clockwise from north). - * - * @param origin Starting point - * @param distance Distance to travel in meters - * @param heading Direction to travel in degrees (0 = North, 90 = East) - * @return New GeoPointInterface position - */ - fun computeOffset( - origin: GeoPointInterface, - distance: Double, - heading: Double, - ): GeoPoint { - val distanceRad = distance / Earth.RADIUS_METERS - val headingRad = heading * DEG_TO_RAD - val lat1Rad = origin.latitude * DEG_TO_RAD - val lng1Rad = origin.longitude * DEG_TO_RAD - - val lat2Rad = - asin( - sin(lat1Rad) * cos(distanceRad) + - cos(lat1Rad) * sin(distanceRad) * cos(headingRad), - ) - - val lng2Rad = - lng1Rad + - atan2( - sin(headingRad) * sin(distanceRad) * cos(lat1Rad), - cos(distanceRad) - sin(lat1Rad) * sin(lat2Rad), - ) - - return GeoPoint( - latitude = lat2Rad * RAD_TO_DEG, - longitude = lng2Rad * RAD_TO_DEG, - altitude = origin.altitude ?: 0.0, - ) - } - - /** - * Returns the location of origin when provided with a GeoPointInterface destination, - * meters travelled and original heading. - * Headings are expressed in degrees clockwise from North. - * - * @param to Destination point - * @param distance Distance travelled in meters - * @param heading Original heading in degrees - * @return Origin GeoPointInterface position, or null if no solution is available - */ - fun computeOffsetOrigin( - to: GeoPointInterface, - distance: Double, - heading: Double, - ): GeoPoint? { - // Calculate the reverse heading - val reverseHeading = (heading + 180) % 360 - - return try { - computeOffset(to, distance, reverseHeading) - } catch (e: Exception) { - null - } - } - - /** - * Returns the length of the given path in meters. - * - * @param path List of GeoPointInterface locations defining the path - * @return Length in meters - */ - fun computeLength(path: List): Double { - if (path.size < 2) return 0.0 - - var length = 0.0 - for (i in 1 until path.size) { - length += computeDistanceBetween(path[i - 1], path[i]) - } - - return length - } - - /** - * Returns the area of a closed path in square meters. - * The path should form a closed polygon. - * - * @param path List of GeoPointInterface locations defining the closed path - * @return Area in square meters - */ - fun computeArea(path: List): Double = abs(computeSignedArea(path)) - - /** - * Returns the signed area of a closed path in square meters. - * The signed area may be used to determine the orientation of the path. - * Positive values indicate counter-clockwise orientation, - * negative values indicate clockwise orientation. - * - * @param path List of GeoPointInterface locations defining the closed path - * @return Signed area in square meters - */ - fun computeSignedArea(path: List): Double { - if (path.size < 3) return 0.0 - - var area = 0.0 - val pointCount = path.size - - for (i in path.indices) { - val j = (i + 1) % pointCount - val lat1 = path[i].latitude * DEG_TO_RAD - val lat2 = path[j].latitude * DEG_TO_RAD - val deltaLng = (path[j].longitude - path[i].longitude) * DEG_TO_RAD - - area += deltaLng * (2 + sin(lat1) + sin(lat2)) - } - - return area * Earth.RADIUS_METERS * Earth.RADIUS_METERS / 2.0 - } - - /** - * Interpolates between two GeoPointInterface locations along the great circle path using spherical linear interpolation (Slerp). - * This method considers Earth's curvature and provides high accuracy for any distance. - * - * @param from Starting point - * @param to Ending point - * @param fraction Interpolation fraction (0.0 = from, 1.0 = to) - * @return Interpolated GeoPointInterface position - */ - fun interpolate( - from: GeoPointInterface, - to: GeoPointInterface, - fraction: Double, - ): GeoPoint { - // ラジアンに変換 - val lat1 = from.latitude * DEG_TO_RAD - val lng1 = from.longitude * DEG_TO_RAD - val lat2 = to.latitude * DEG_TO_RAD - val lng2 = to.longitude * DEG_TO_RAD - - // 3D単位ベクトルに変換 - val x1 = cos(lat1) * cos(lng1) - val y1 = cos(lat1) * sin(lng1) - val z1 = sin(lat1) - - val x2 = cos(lat2) * cos(lng2) - val y2 = cos(lat2) * sin(lng2) - val z2 = sin(lat2) - - // 内積から角度を求める - val dot = x1 * x2 + y1 * y2 + z1 * z2 - val angle = acos(dot.coerceIn(-1.0, 1.0)) - - // 非常に近い点は線形補間 - if (angle < 1e-6) { - val interpolatedAltitude = - when { - from.altitude != null && to.altitude != null -> - from.altitude!! + fraction * (to.altitude!! - from.altitude!!) - from.altitude != null -> from.altitude - to.altitude != null -> to.altitude - else -> 0.0 - } - - return GeoPoint( - latitude = from.latitude + fraction * (to.latitude - from.latitude), - longitude = from.longitude + fraction * (to.longitude - from.longitude), - altitude = interpolatedAltitude!!, - ) - } - - // 球面線形補間（Slerp） - val sinAngle = sin(angle) - val weightFrom = sin((1 - fraction) * angle) / sinAngle - val weightTo = sin(fraction * angle) / sinAngle - - val vectorX = weightFrom * x1 + weightTo * x2 - val vectorY = weightFrom * y1 + weightTo * y2 - val vectorZ = weightFrom * z1 + weightTo * z2 - - // 3Dベクトルから緯度経度に変換 - val lat = asin(vectorZ) * RAD_TO_DEG - val lng = atan2(vectorY, vectorX) * RAD_TO_DEG - - val interpolatedAltitude = - when { - from.altitude != null && to.altitude != null -> - from.altitude!! + fraction * (to.altitude!! - from.altitude!!) - from.altitude != null -> from.altitude - to.altitude != null -> to.altitude - else -> 0.0 - } - - return GeoPoint( - latitude = lat, - longitude = lng, - altitude = interpolatedAltitude!!, - ) - } - - // `linearInterpolate` moved to `Planar.interpolate` (the straight-line model)。 - // 説明もそちらの KDoc にある。 - - /** - * Clamps latitude to the range [-90, 90]. - */ - private fun clampLat(lat: Double): Double = lat.coerceIn(-90.0, 90.0) -} +package com.mapconductor.core.spherical + +import com.mapconductor.core.features.GeoPoint +import com.mapconductor.core.features.GeoPointInterface +import com.mapconductor.core.projection.Earth +import kotlin.math.abs +import kotlin.math.acos +import kotlin.math.asin +import kotlin.math.atan2 +import kotlin.math.cos +import kotlin.math.sin +import kotlin.math.sqrt + +/** + * Spherical geometry utility functions for calculating distances, headings, + * and positions on Earth's surface using the spherical Earth model. + * + * This is a Kotlin port of the Cordova Google Maps plugin spherical.js utilities. + * Uses GeoPointInterface instead of LatLng for coordinate representation. + */ +object Spherical { + // Mathematical constants + private const val PI = Math.PI + private const val RAD_TO_DEG = 180.0 / PI + private const val DEG_TO_RAD = PI / 180.0 + + /** + * Returns the distance, in meters, between two GeoPointInterface locations. + * Uses the haversine formula. + * + * @param from Starting point + * @param to Ending point + * @return Distance in meters + */ + fun computeDistanceBetween( + from: GeoPointInterface, + to: GeoPointInterface, + ): Double { + val lat1Rad = from.latitude * DEG_TO_RAD + val lat2Rad = to.latitude * DEG_TO_RAD + val deltaLat = (to.latitude - from.latitude) * DEG_TO_RAD + val deltaLng = (to.longitude - from.longitude) * DEG_TO_RAD + + val haversineA = + sin(deltaLat / 2) * sin(deltaLat / 2) + + cos(lat1Rad) * cos(lat2Rad) * + sin(deltaLng / 2) * sin(deltaLng / 2) + + val centralAngle = 2 * atan2(sqrt(haversineA), sqrt(1 - haversineA)) + + return Earth.RADIUS_METERS * centralAngle + } + + /** + * Returns the heading from one GeoPointInterface to another GeoPointInterface. + * Headings are expressed in degrees clockwise from North within the range (-180, 180). + * + * @param from Starting point + * @param to Ending point + * @return Heading in degrees + */ + fun computeHeading( + from: GeoPointInterface, + to: GeoPointInterface, + ): Double { + val lat1Rad = from.latitude * DEG_TO_RAD + val lat2Rad = to.latitude * DEG_TO_RAD + val deltaLng = (to.longitude - from.longitude) * DEG_TO_RAD + + val deltaY = sin(deltaLng) * cos(lat2Rad) + val deltaX = cos(lat1Rad) * sin(lat2Rad) - sin(lat1Rad) * cos(lat2Rad) * cos(deltaLng) + + var heading = atan2(deltaY, deltaX) * RAD_TO_DEG + + // Normalize to (-180, 180] + while (heading > 180) heading -= 360 + while (heading <= -180) heading += 360 + + return heading + } + + /** + * Returns the GeoPointInterface resulting from moving a distance from an origin + * in the specified heading (expressed in degrees clockwise from north). + * + * @param origin Starting point + * @param distance Distance to travel in meters + * @param heading Direction to travel in degrees (0 = North, 90 = East) + * @return New GeoPointInterface position + */ + fun computeOffset( + origin: GeoPointInterface, + distance: Double, + heading: Double, + ): GeoPoint { + val distanceRad = distance / Earth.RADIUS_METERS + val headingRad = heading * DEG_TO_RAD + val lat1Rad = origin.latitude * DEG_TO_RAD + val lng1Rad = origin.longitude * DEG_TO_RAD + + val lat2Rad = + asin( + sin(lat1Rad) * cos(distanceRad) + + cos(lat1Rad) * sin(distanceRad) * cos(headingRad), + ) + + val lng2Rad = + lng1Rad + + atan2( + sin(headingRad) * sin(distanceRad) * cos(lat1Rad), + cos(distanceRad) - sin(lat1Rad) * sin(lat2Rad), + ) + + return GeoPoint( + latitude = lat2Rad * RAD_TO_DEG, + longitude = lng2Rad * RAD_TO_DEG, + altitude = origin.altitude ?: 0.0, + ) + } + + /** + * Returns the location of origin when provided with a GeoPointInterface destination, + * meters travelled and original heading. + * Headings are expressed in degrees clockwise from North. + * + * @param to Destination point + * @param distance Distance travelled in meters + * @param heading Original heading in degrees + * @return Origin GeoPointInterface position, or null if no solution is available + */ + fun computeOffsetOrigin( + to: GeoPointInterface, + distance: Double, + heading: Double, + ): GeoPoint? { + // Calculate the reverse heading + val reverseHeading = (heading + 180) % 360 + + return try { + computeOffset(to, distance, reverseHeading) + } catch (e: Exception) { + null + } + } + + /** + * Returns the length of the given path in meters. + * + * @param path List of GeoPointInterface locations defining the path + * @return Length in meters + */ + fun computeLength(path: List): Double { + if (path.size < 2) return 0.0 + + var length = 0.0 + for (i in 1 until path.size) { + length += computeDistanceBetween(path[i - 1], path[i]) + } + + return length + } + + /** + * Returns the area of a closed path in square meters. + * The path should form a closed polygon. + * + * @param path List of GeoPointInterface locations defining the closed path + * @return Area in square meters + */ + fun computeArea(path: List): Double = abs(computeSignedArea(path)) + + /** + * Returns the signed area of a closed path in square meters. + * The signed area may be used to determine the orientation of the path. + * Positive values indicate clockwise orientation, + * negative values indicate counter-clockwise orientation + * (as seen on a north-up map, e.g. SW -> NW -> NE -> SE is clockwise and positive). + * + * @param path List of GeoPointInterface locations defining the closed path + * @return Signed area in square meters + */ + fun computeSignedArea(path: List): Double { + if (path.size < 3) return 0.0 + + var area = 0.0 + val pointCount = path.size + + for (i in path.indices) { + val j = (i + 1) % pointCount + val lat1 = path[i].latitude * DEG_TO_RAD + val lat2 = path[j].latitude * DEG_TO_RAD + val deltaLng = (path[j].longitude - path[i].longitude) * DEG_TO_RAD + + area += deltaLng * (2 + sin(lat1) + sin(lat2)) + } + + return area * Earth.RADIUS_METERS * Earth.RADIUS_METERS / 2.0 + } + + /** + * Interpolates between two GeoPointInterface locations along the great circle path using spherical linear interpolation (Slerp). + * This method considers Earth's curvature and provides high accuracy for any distance. + * + * @param from Starting point + * @param to Ending point + * @param fraction Interpolation fraction (0.0 = from, 1.0 = to) + * @return Interpolated GeoPointInterface position + */ + fun interpolate( + from: GeoPointInterface, + to: GeoPointInterface, + fraction: Double, + ): GeoPoint { + // ラジアンに変換 + val lat1 = from.latitude * DEG_TO_RAD + val lng1 = from.longitude * DEG_TO_RAD + val lat2 = to.latitude * DEG_TO_RAD + val lng2 = to.longitude * DEG_TO_RAD + + // 3D単位ベクトルに変換 + val x1 = cos(lat1) * cos(lng1) + val y1 = cos(lat1) * sin(lng1) + val z1 = sin(lat1) + + val x2 = cos(lat2) * cos(lng2) + val y2 = cos(lat2) * sin(lng2) + val z2 = sin(lat2) + + // 内積から角度を求める + val dot = x1 * x2 + y1 * y2 + z1 * z2 + val angle = acos(dot.coerceIn(-1.0, 1.0)) + + // 非常に近い点は線形補間 + if (angle < 1e-6) { + val interpolatedAltitude = + when { + from.altitude != null && to.altitude != null -> + from.altitude!! + fraction * (to.altitude!! - from.altitude!!) + from.altitude != null -> from.altitude + to.altitude != null -> to.altitude + else -> 0.0 + } + + return GeoPoint( + latitude = from.latitude + fraction * (to.latitude - from.latitude), + longitude = from.longitude + fraction * (to.longitude - from.longitude), + altitude = interpolatedAltitude!!, + ) + } + + // 球面線形補間（Slerp） + val sinAngle = sin(angle) + val weightFrom = sin((1 - fraction) * angle) / sinAngle + val weightTo = sin(fraction * angle) / sinAngle + + val vectorX = weightFrom * x1 + weightTo * x2 + val vectorY = weightFrom * y1 + weightTo * y2 + val vectorZ = weightFrom * z1 + weightTo * z2 + + // 3Dベクトルから緯度経度に変換 + val lat = asin(vectorZ) * RAD_TO_DEG + val lng = atan2(vectorY, vectorX) * RAD_TO_DEG + + val interpolatedAltitude = + when { + from.altitude != null && to.altitude != null -> + from.altitude!! + fraction * (to.altitude!! - from.altitude!!) + from.altitude != null -> from.altitude + to.altitude != null -> to.altitude + else -> 0.0 + } + + return GeoPoint( + latitude = lat, + longitude = lng, + altitude = interpolatedAltitude!!, + ) + } + + // `linearInterpolate` moved to `Planar.interpolate` (the straight-line model)。 + // 説明もそちらの KDoc にある。 + + /** + * Clamps latitude to the range [-90, 90]. + */ + private fun clampLat(lat: Double): Double = lat.coerceIn(-90.0, 90.0) +} diff --git a/src/main/java/com/mapconductor/core/spherical/WGS84Geodesic.kt b/src/main/java/com/mapconductor/core/spherical/WGS84Geodesic.kt index 35a3521..198dd01 100644 --- a/src/main/java/com/mapconductor/core/spherical/WGS84Geodesic.kt +++ b/src/main/java/com/mapconductor/core/spherical/WGS84Geodesic.kt @@ -1,323 +1,326 @@ -package com.mapconductor.core.spherical - -import com.mapconductor.core.features.GeoPoint -import com.mapconductor.core.features.GeoPointInterface -import com.mapconductor.core.projection.Earth -import kotlin.math.abs -import kotlin.math.atan -import kotlin.math.atan2 -import kotlin.math.cos -import kotlin.math.ln -import kotlin.math.sin -import kotlin.math.sqrt - -/** - * WGS84 ellipsoid geodesic distance/interpolation via Vincenty's formulae. - * Ported from the react-sdk/ios-sdk implementations so all platforms agree - * without depending on an external geographiclib library. Vincenty can fail - * to converge for near-antipodal point pairs; in that (rare, map-rendering- - * irrelevant) case this falls back to a spherical (haversine) approximation, - * matching react/ios behavior exactly. - */ -object WGS84Geodesic { - private const val FLATTENING = Earth.FLATTENING - private const val SEMI_MAJOR_AXIS = Earth.RADIUS_METERS - private const val SEMI_MINOR_AXIS = Earth.SEMI_MINOR_AXIS_METERS - - // Authalic-sphere constants for ellipsoidal polygon area (see computeSignedArea). - private val ECCENTRICITY_SQ = FLATTENING * (2 - FLATTENING) - private val ECCENTRICITY = sqrt(ECCENTRICITY_SQ) - private val AUTHALIC_QP = - 1 - ((1 - ECCENTRICITY_SQ) / (2 * ECCENTRICITY)) * - ln((1 - ECCENTRICITY) / (1 + ECCENTRICITY)) - private val AUTHALIC_RADIUS = SEMI_MAJOR_AXIS * sqrt(AUTHALIC_QP / 2) - - private data class InverseResult( - val distanceMeters: Double, - val initialBearingRad: Double, - ) - - fun computeDistanceBetween( - from: GeoPointInterface, - to: GeoPointInterface, - ): Double = inverse(from, to).distanceMeters - - fun computeHeading( - from: GeoPointInterface, - to: GeoPointInterface, - ): Double = Math.toDegrees(inverse(from, to).initialBearingRad) - - fun computeOffset( - origin: GeoPointInterface, - distance: Double, - heading: Double, - ): GeoPoint { - val destination = direct(origin, Math.toRadians(heading), distance) - return GeoPoint(destination.first, destination.second, origin.altitude ?: 0.0) - } - - fun computeOffsetOrigin( - to: GeoPointInterface, - distance: Double, - heading: Double, - ): GeoPoint? = computeOffset(to, distance, (heading + 180.0) % 360.0) - - fun computeLength(path: List): Double { - if (path.size < 2) return 0.0 - var length = 0.0 - for (i in 1 until path.size) { - length += computeDistanceBetween(path[i - 1], path[i]) - } - return length - } - - /** Densify a path by inserting points along the geodesic. */ - fun createInterpolatePoints( - points: List, - maxSegmentLength: Double = 10000.0, - ): List = densifyAlongGeodesic(points, maxSegmentLength) - - /** Closest point on the geodesic segment within [thresholdMeters], or null. */ - fun pointOnLineOrNull( - from: GeoPointInterface, - to: GeoPointInterface, - position: GeoPointInterface, - thresholdMeters: Double, - ): Pair? = geodesicPointOnLineOrNull(from, to, position, thresholdMeters) - - /** - * Ellipsoidal (authalic-sphere) polygon area. Maps geodetic latitudes to - * authalic latitudes on the equal-area sphere of radius [AUTHALIC_RADIUS] - * (which reproduces the exact WGS84 surface area), then applies the - * spherical-excess formula. Accounts for flattening — unlike the sphere-based - * [Spherical.computeSignedArea] — matching the ellipsoidal area to well under - * 0.01% for typical polygons. Edges are treated as authalic-sphere arcs; the - * higher-order geodesic-edge terms of Karney's exact method are omitted. - */ - fun computeSignedArea(path: List): Double { - if (path.size < 3) return 0.0 - var area = 0.0 - val pointCount = path.size - for (i in 0 until pointCount) { - val j = (i + 1) % pointCount - val sinXi1 = authalicSinLatitude(path[i].latitude) - val sinXi2 = authalicSinLatitude(path[j].latitude) - val deltaLng = Math.toRadians(path[j].longitude - path[i].longitude) - area += deltaLng * (2 + sinXi1 + sinXi2) - } - return area * AUTHALIC_RADIUS * AUTHALIC_RADIUS / 2.0 - } - - fun computeArea(path: List): Double = abs(computeSignedArea(path)) - - private fun authalicSinLatitude(latitudeDeg: Double): Double { - val sinPhi = sin(Math.toRadians(latitudeDeg)) - val q = - (1 - ECCENTRICITY_SQ) * - ( - sinPhi / (1 - ECCENTRICITY_SQ * sinPhi * sinPhi) - - (1 / (2 * ECCENTRICITY)) * - ln((1 - ECCENTRICITY * sinPhi) / (1 + ECCENTRICITY * sinPhi)) - ) - return q / AUTHALIC_QP - } - - fun interpolate( - from: GeoPointInterface, - to: GeoPointInterface, - fraction: Double, - ): GeoPoint { - val line = inverse(from, to) - val destination = direct(from, line.initialBearingRad, line.distanceMeters * fraction) - - val altitude = - when { - from.altitude != null && to.altitude != null -> - from.altitude!! + fraction * (to.altitude!! - from.altitude!!) - from.altitude != null -> from.altitude!! - to.altitude != null -> to.altitude!! - else -> 0.0 - } - - return GeoPoint(destination.first, destination.second, altitude) - } - - private fun inverse( - from: GeoPointInterface, - to: GeoPointInterface, - ): InverseResult { - val lat1 = Math.toRadians(from.latitude) - val lat2 = Math.toRadians(to.latitude) - val lon1 = Math.toRadians(from.longitude) - val lon2 = Math.toRadians(to.longitude) - val longitudeDifference = lon2 - lon1 - - val u1 = atan((1 - FLATTENING) * Math.tan(lat1)) - val u2 = atan((1 - FLATTENING) * Math.tan(lat2)) - val sinU1 = sin(u1) - val cosU1 = cos(u1) - val sinU2 = sin(u2) - val cosU2 = cos(u2) - - var lambda = longitudeDifference - var lambdaPrev: Double - var iterLimit = 100 - var cosSqAlpha = 0.0 - var sinSigma = 0.0 - var cos2SigmaM = 0.0 - var cosSigma = 0.0 - var sigma = 0.0 - - do { - val sinLambda = sin(lambda) - val cosLambda = cos(lambda) - sinSigma = - sqrt( - (cosU2 * sinLambda) * (cosU2 * sinLambda) + - (cosU1 * sinU2 - sinU1 * cosU2 * cosLambda) * - (cosU1 * sinU2 - sinU1 * cosU2 * cosLambda), - ) - - if (sinSigma == 0.0) { - return InverseResult(0.0, 0.0) - } - - cosSigma = sinU1 * sinU2 + cosU1 * cosU2 * cosLambda - sigma = atan2(sinSigma, cosSigma) - val sinAlpha = cosU1 * cosU2 * sinLambda / sinSigma - cosSqAlpha = 1 - sinAlpha * sinAlpha - cos2SigmaM = cosSigma - 2 * sinU1 * sinU2 / cosSqAlpha - if (!cos2SigmaM.isFinite()) cos2SigmaM = 0.0 - - val correctionFactor = FLATTENING / 16 * cosSqAlpha * (4 + FLATTENING * (4 - 3 * cosSqAlpha)) - lambdaPrev = lambda - lambda = - longitudeDifference + - (1 - correctionFactor) * FLATTENING * sinAlpha * - ( - sigma + - correctionFactor * sinSigma * - (cos2SigmaM + correctionFactor * cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM)) - ) - iterLimit -= 1 - } while (abs(lambda - lambdaPrev) > 1e-12 && iterLimit > 0) - - if (iterLimit == 0) return sphericalFallbackInverse(from, to) - - val uSq = - cosSqAlpha * (SEMI_MAJOR_AXIS * SEMI_MAJOR_AXIS - SEMI_MINOR_AXIS * SEMI_MINOR_AXIS) / - (SEMI_MINOR_AXIS * SEMI_MINOR_AXIS) - val ellipsoidFactor = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq))) - val correctionTerm = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq))) - val deltaSigma = - correctionTerm * sinSigma * ( - cos2SigmaM + correctionTerm / 4 * ( - cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM) - - correctionTerm / 6 * cos2SigmaM * (-3 + 4 * sinSigma * sinSigma) * - (-3 + 4 * cos2SigmaM * cos2SigmaM) - ) - ) - - val distance = SEMI_MINOR_AXIS * ellipsoidFactor * (sigma - deltaSigma) - val initialBearing = - atan2( - cosU2 * sin(lambda), - cosU1 * sinU2 - sinU1 * cosU2 * cos(lambda), - ) - - return InverseResult(distance, initialBearing) - } - - private fun direct( - origin: GeoPointInterface, - initialBearingRad: Double, - distanceMeters: Double, - ): Pair { - val lat1 = Math.toRadians(origin.latitude) - val lon1 = Math.toRadians(origin.longitude) - val sinAlpha1 = sin(initialBearingRad) - val cosAlpha1 = cos(initialBearingRad) - - val tanU1 = (1 - FLATTENING) * Math.tan(lat1) - val cosU1 = 1 / sqrt(1 + tanU1 * tanU1) - val sinU1 = tanU1 * cosU1 - val sigma1 = atan2(tanU1, cosAlpha1) - val sinAlpha = cosU1 * sinAlpha1 - val cosSqAlpha = 1 - sinAlpha * sinAlpha - val uSq = - cosSqAlpha * (SEMI_MAJOR_AXIS * SEMI_MAJOR_AXIS - SEMI_MINOR_AXIS * SEMI_MINOR_AXIS) / - (SEMI_MINOR_AXIS * SEMI_MINOR_AXIS) - val ellipsoidFactor = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq))) - val correctionTerm = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq))) - - var sigma = distanceMeters / (SEMI_MINOR_AXIS * ellipsoidFactor) - var sigmaPrev: Double - var cos2SigmaM = 0.0 - var sinSigma = 0.0 - var cosSigma = 0.0 - // ios-sdk と同じく非収束の病的入力での無限ループを防ぐ反復上限。 - var iterLimit = 1000 - - do { - cos2SigmaM = cos(2 * sigma1 + sigma) - sinSigma = sin(sigma) - cosSigma = cos(sigma) - val deltaSigma = - correctionTerm * sinSigma * ( - cos2SigmaM + correctionTerm / 4 * ( - cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM) - - correctionTerm / 6 * cos2SigmaM * (-3 + 4 * sinSigma * sinSigma) * - (-3 + 4 * cos2SigmaM * cos2SigmaM) - ) - ) - sigmaPrev = sigma - sigma = distanceMeters / (SEMI_MINOR_AXIS * ellipsoidFactor) + deltaSigma - iterLimit-- - } while (abs(sigma - sigmaPrev) > 1e-12 && iterLimit > 0) - - val tmp = sinU1 * sinSigma - cosU1 * cosSigma * cosAlpha1 - val lat2 = - atan2( - sinU1 * cosSigma + cosU1 * sinSigma * cosAlpha1, - (1 - FLATTENING) * sqrt(sinAlpha * sinAlpha + tmp * tmp), - ) - val lambda = - atan2( - sinSigma * sinAlpha1, - cosU1 * cosSigma - sinU1 * sinSigma * cosAlpha1, - ) - val correctionFactor = FLATTENING / 16 * cosSqAlpha * (4 + FLATTENING * (4 - 3 * cosSqAlpha)) - val longitudeDifference = - lambda - (1 - correctionFactor) * FLATTENING * sinAlpha * - ( - sigma + - correctionFactor * sinSigma * - (cos2SigmaM + correctionFactor * cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM)) - ) - val lon2 = lon1 + longitudeDifference - - return Pair(Math.toDegrees(lat2), normalizeLng(Math.toDegrees(lon2))) - } - - private fun sphericalFallbackInverse( - from: GeoPointInterface, - to: GeoPointInterface, - ): InverseResult { - val lat1 = Math.toRadians(from.latitude) - val lat2 = Math.toRadians(to.latitude) - val deltaLat = Math.toRadians(to.latitude - from.latitude) - val deltaLng = Math.toRadians(to.longitude - from.longitude) - val a = - sin(deltaLat / 2) * sin(deltaLat / 2) + - cos(lat1) * cos(lat2) * sin(deltaLng / 2) * sin(deltaLng / 2) - val centralAngle = 2 * atan2(sqrt(a), sqrt(1 - a)) - val y = sin(deltaLng) * cos(lat2) - val x = cos(lat1) * sin(lat2) - sin(lat1) * cos(lat2) * cos(deltaLng) - return InverseResult( - distanceMeters = SEMI_MAJOR_AXIS * centralAngle, - initialBearingRad = atan2(y, x), - ) - } - - private fun normalizeLng(lng: Double): Double = ((((lng + 180.0) % 360.0) + 360.0) % 360.0) - 180.0 -} +package com.mapconductor.core.spherical + +import com.mapconductor.core.features.GeoPoint +import com.mapconductor.core.features.GeoPointInterface +import com.mapconductor.core.projection.Earth +import kotlin.math.abs +import kotlin.math.atan +import kotlin.math.atan2 +import kotlin.math.cos +import kotlin.math.ln +import kotlin.math.sin +import kotlin.math.sqrt + +/** + * WGS84 ellipsoid geodesic distance/interpolation via Vincenty's formulae. + * Ported from the react-sdk/ios-sdk implementations so all platforms agree + * without depending on an external geographiclib library. Vincenty can fail + * to converge for near-antipodal point pairs; in that (rare, map-rendering- + * irrelevant) case this falls back to a spherical (haversine) approximation, + * matching react/ios behavior exactly. + */ +object WGS84Geodesic { + private const val FLATTENING = Earth.FLATTENING + private const val SEMI_MAJOR_AXIS = Earth.RADIUS_METERS + private const val SEMI_MINOR_AXIS = Earth.SEMI_MINOR_AXIS_METERS + + // Authalic-sphere constants for ellipsoidal polygon area (see computeSignedArea). + private val ECCENTRICITY_SQ = FLATTENING * (2 - FLATTENING) + private val ECCENTRICITY = sqrt(ECCENTRICITY_SQ) + private val AUTHALIC_QP = + 1 - ((1 - ECCENTRICITY_SQ) / (2 * ECCENTRICITY)) * + ln((1 - ECCENTRICITY) / (1 + ECCENTRICITY)) + private val AUTHALIC_RADIUS = SEMI_MAJOR_AXIS * sqrt(AUTHALIC_QP / 2) + + private data class InverseResult( + val distanceMeters: Double, + val initialBearingRad: Double, + ) + + fun computeDistanceBetween( + from: GeoPointInterface, + to: GeoPointInterface, + ): Double = inverse(from, to).distanceMeters + + fun computeHeading( + from: GeoPointInterface, + to: GeoPointInterface, + ): Double = Math.toDegrees(inverse(from, to).initialBearingRad) + + fun computeOffset( + origin: GeoPointInterface, + distance: Double, + heading: Double, + ): GeoPoint { + val destination = direct(origin, Math.toRadians(heading), distance) + return GeoPoint(destination.first, destination.second, origin.altitude ?: 0.0) + } + + fun computeOffsetOrigin( + to: GeoPointInterface, + distance: Double, + heading: Double, + ): GeoPoint? = computeOffset(to, distance, (heading + 180.0) % 360.0) + + fun computeLength(path: List): Double { + if (path.size < 2) return 0.0 + var length = 0.0 + for (i in 1 until path.size) { + length += computeDistanceBetween(path[i - 1], path[i]) + } + return length + } + + /** Densify a path by inserting points along the geodesic. */ + fun createInterpolatePoints( + points: List, + maxSegmentLength: Double = 10000.0, + ): List = densifyAlongGeodesic(points, maxSegmentLength) + + /** Closest point on the geodesic segment within [thresholdMeters], or null. */ + fun pointOnLineOrNull( + from: GeoPointInterface, + to: GeoPointInterface, + position: GeoPointInterface, + thresholdMeters: Double, + ): Pair? = geodesicPointOnLineOrNull(from, to, position, thresholdMeters) + + /** + * Ellipsoidal (authalic-sphere) polygon area. Maps geodetic latitudes to + * authalic latitudes on the equal-area sphere of radius [AUTHALIC_RADIUS] + * (which reproduces the exact WGS84 surface area), then applies the + * spherical-excess formula. Accounts for flattening — unlike the sphere-based + * [Spherical.computeSignedArea] — matching the ellipsoidal area to well under + * 0.01% for typical polygons. Edges are treated as authalic-sphere arcs; the + * higher-order geodesic-edge terms of Karney's exact method are omitted. + * + * The sign follows [Spherical.computeSignedArea]: positive for a clockwise + * ring, negative for a counter-clockwise ring (as seen on a north-up map). + */ + fun computeSignedArea(path: List): Double { + if (path.size < 3) return 0.0 + var area = 0.0 + val pointCount = path.size + for (i in 0 until pointCount) { + val j = (i + 1) % pointCount + val sinXi1 = authalicSinLatitude(path[i].latitude) + val sinXi2 = authalicSinLatitude(path[j].latitude) + val deltaLng = Math.toRadians(path[j].longitude - path[i].longitude) + area += deltaLng * (2 + sinXi1 + sinXi2) + } + return area * AUTHALIC_RADIUS * AUTHALIC_RADIUS / 2.0 + } + + fun computeArea(path: List): Double = abs(computeSignedArea(path)) + + private fun authalicSinLatitude(latitudeDeg: Double): Double { + val sinPhi = sin(Math.toRadians(latitudeDeg)) + val q = + (1 - ECCENTRICITY_SQ) * + ( + sinPhi / (1 - ECCENTRICITY_SQ * sinPhi * sinPhi) - + (1 / (2 * ECCENTRICITY)) * + ln((1 - ECCENTRICITY * sinPhi) / (1 + ECCENTRICITY * sinPhi)) + ) + return q / AUTHALIC_QP + } + + fun interpolate( + from: GeoPointInterface, + to: GeoPointInterface, + fraction: Double, + ): GeoPoint { + val line = inverse(from, to) + val destination = direct(from, line.initialBearingRad, line.distanceMeters * fraction) + + val altitude = + when { + from.altitude != null && to.altitude != null -> + from.altitude!! + fraction * (to.altitude!! - from.altitude!!) + from.altitude != null -> from.altitude!! + to.altitude != null -> to.altitude!! + else -> 0.0 + } + + return GeoPoint(destination.first, destination.second, altitude) + } + + private fun inverse( + from: GeoPointInterface, + to: GeoPointInterface, + ): InverseResult { + val lat1 = Math.toRadians(from.latitude) + val lat2 = Math.toRadians(to.latitude) + val lon1 = Math.toRadians(from.longitude) + val lon2 = Math.toRadians(to.longitude) + val longitudeDifference = lon2 - lon1 + + val u1 = atan((1 - FLATTENING) * Math.tan(lat1)) + val u2 = atan((1 - FLATTENING) * Math.tan(lat2)) + val sinU1 = sin(u1) + val cosU1 = cos(u1) + val sinU2 = sin(u2) + val cosU2 = cos(u2) + + var lambda = longitudeDifference + var lambdaPrev: Double + var iterLimit = 100 + var cosSqAlpha = 0.0 + var sinSigma = 0.0 + var cos2SigmaM = 0.0 + var cosSigma = 0.0 + var sigma = 0.0 + + do { + val sinLambda = sin(lambda) + val cosLambda = cos(lambda) + sinSigma = + sqrt( + (cosU2 * sinLambda) * (cosU2 * sinLambda) + + (cosU1 * sinU2 - sinU1 * cosU2 * cosLambda) * + (cosU1 * sinU2 - sinU1 * cosU2 * cosLambda), + ) + + if (sinSigma == 0.0) { + return InverseResult(0.0, 0.0) + } + + cosSigma = sinU1 * sinU2 + cosU1 * cosU2 * cosLambda + sigma = atan2(sinSigma, cosSigma) + val sinAlpha = cosU1 * cosU2 * sinLambda / sinSigma + cosSqAlpha = 1 - sinAlpha * sinAlpha + cos2SigmaM = cosSigma - 2 * sinU1 * sinU2 / cosSqAlpha + if (!cos2SigmaM.isFinite()) cos2SigmaM = 0.0 + + val correctionFactor = FLATTENING / 16 * cosSqAlpha * (4 + FLATTENING * (4 - 3 * cosSqAlpha)) + lambdaPrev = lambda + lambda = + longitudeDifference + + (1 - correctionFactor) * FLATTENING * sinAlpha * + ( + sigma + + correctionFactor * sinSigma * + (cos2SigmaM + correctionFactor * cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM)) + ) + iterLimit -= 1 + } while (abs(lambda - lambdaPrev) > 1e-12 && iterLimit > 0) + + if (iterLimit == 0) return sphericalFallbackInverse(from, to) + + val uSq = + cosSqAlpha * (SEMI_MAJOR_AXIS * SEMI_MAJOR_AXIS - SEMI_MINOR_AXIS * SEMI_MINOR_AXIS) / + (SEMI_MINOR_AXIS * SEMI_MINOR_AXIS) + val ellipsoidFactor = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq))) + val correctionTerm = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq))) + val deltaSigma = + correctionTerm * sinSigma * ( + cos2SigmaM + correctionTerm / 4 * ( + cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM) - + correctionTerm / 6 * cos2SigmaM * (-3 + 4 * sinSigma * sinSigma) * + (-3 + 4 * cos2SigmaM * cos2SigmaM) + ) + ) + + val distance = SEMI_MINOR_AXIS * ellipsoidFactor * (sigma - deltaSigma) + val initialBearing = + atan2( + cosU2 * sin(lambda), + cosU1 * sinU2 - sinU1 * cosU2 * cos(lambda), + ) + + return InverseResult(distance, initialBearing) + } + + private fun direct( + origin: GeoPointInterface, + initialBearingRad: Double, + distanceMeters: Double, + ): Pair { + val lat1 = Math.toRadians(origin.latitude) + val lon1 = Math.toRadians(origin.longitude) + val sinAlpha1 = sin(initialBearingRad) + val cosAlpha1 = cos(initialBearingRad) + + val tanU1 = (1 - FLATTENING) * Math.tan(lat1) + val cosU1 = 1 / sqrt(1 + tanU1 * tanU1) + val sinU1 = tanU1 * cosU1 + val sigma1 = atan2(tanU1, cosAlpha1) + val sinAlpha = cosU1 * sinAlpha1 + val cosSqAlpha = 1 - sinAlpha * sinAlpha + val uSq = + cosSqAlpha * (SEMI_MAJOR_AXIS * SEMI_MAJOR_AXIS - SEMI_MINOR_AXIS * SEMI_MINOR_AXIS) / + (SEMI_MINOR_AXIS * SEMI_MINOR_AXIS) + val ellipsoidFactor = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq))) + val correctionTerm = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq))) + + var sigma = distanceMeters / (SEMI_MINOR_AXIS * ellipsoidFactor) + var sigmaPrev: Double + var cos2SigmaM = 0.0 + var sinSigma = 0.0 + var cosSigma = 0.0 + // ios-sdk と同じく非収束の病的入力での無限ループを防ぐ反復上限。 + var iterLimit = 1000 + + do { + cos2SigmaM = cos(2 * sigma1 + sigma) + sinSigma = sin(sigma) + cosSigma = cos(sigma) + val deltaSigma = + correctionTerm * sinSigma * ( + cos2SigmaM + correctionTerm / 4 * ( + cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM) - + correctionTerm / 6 * cos2SigmaM * (-3 + 4 * sinSigma * sinSigma) * + (-3 + 4 * cos2SigmaM * cos2SigmaM) + ) + ) + sigmaPrev = sigma + sigma = distanceMeters / (SEMI_MINOR_AXIS * ellipsoidFactor) + deltaSigma + iterLimit-- + } while (abs(sigma - sigmaPrev) > 1e-12 && iterLimit > 0) + + val tmp = sinU1 * sinSigma - cosU1 * cosSigma * cosAlpha1 + val lat2 = + atan2( + sinU1 * cosSigma + cosU1 * sinSigma * cosAlpha1, + (1 - FLATTENING) * sqrt(sinAlpha * sinAlpha + tmp * tmp), + ) + val lambda = + atan2( + sinSigma * sinAlpha1, + cosU1 * cosSigma - sinU1 * sinSigma * cosAlpha1, + ) + val correctionFactor = FLATTENING / 16 * cosSqAlpha * (4 + FLATTENING * (4 - 3 * cosSqAlpha)) + val longitudeDifference = + lambda - (1 - correctionFactor) * FLATTENING * sinAlpha * + ( + sigma + + correctionFactor * sinSigma * + (cos2SigmaM + correctionFactor * cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM)) + ) + val lon2 = lon1 + longitudeDifference + + return Pair(Math.toDegrees(lat2), normalizeLng(Math.toDegrees(lon2))) + } + + private fun sphericalFallbackInverse( + from: GeoPointInterface, + to: GeoPointInterface, + ): InverseResult { + val lat1 = Math.toRadians(from.latitude) + val lat2 = Math.toRadians(to.latitude) + val deltaLat = Math.toRadians(to.latitude - from.latitude) + val deltaLng = Math.toRadians(to.longitude - from.longitude) + val a = + sin(deltaLat / 2) * sin(deltaLat / 2) + + cos(lat1) * cos(lat2) * sin(deltaLng / 2) * sin(deltaLng / 2) + val centralAngle = 2 * atan2(sqrt(a), sqrt(1 - a)) + val y = sin(deltaLng) * cos(lat2) + val x = cos(lat1) * sin(lat2) - sin(lat1) * cos(lat2) * cos(deltaLng) + return InverseResult( + distanceMeters = SEMI_MAJOR_AXIS * centralAngle, + initialBearingRad = atan2(y, x), + ) + } + + private fun normalizeLng(lng: Double): Double = ((((lng + 180.0) % 360.0) + 360.0) % 360.0) - 180.0 +} diff --git a/src/test/java/com/mapconductor/core/ShapeDefaultsTest.kt b/src/test/java/com/mapconductor/core/ShapeDefaultsTest.kt new file mode 100644 index 0000000..5965017 --- /dev/null +++ b/src/test/java/com/mapconductor/core/ShapeDefaultsTest.kt @@ -0,0 +1,60 @@ +package com.mapconductor.core + +import androidx.compose.ui.unit.dp +import com.mapconductor.core.circle.CircleState +import com.mapconductor.core.features.GeoPoint +import com.mapconductor.core.features.GeoRectBounds +import com.mapconductor.core.groundimage.GroundImageState +import com.mapconductor.core.polygon.PolygonState +import com.mapconductor.core.polyline.PolylineState +import org.junit.Assert.assertEquals +import org.junit.Test + +/** + * 図形・地表画像の既定値。コンポーネント形（android-sdk-compose の `Circle(center, ...)` など）と + * state 形で同じ値になるよう 3 プラットフォームで揃えている: + * Circle / Polygon / Polyline の strokeWidth = 1.dp、GroundImage の opacity = 1.0f。 + * js-sdk-core の `ShapeDefaults.test.mjs`、ios-sdk の `ShapeDefaultsTests` と同じ観点。 + */ +class ShapeDefaultsTest { + private val p1 = GeoPoint.fromLatLong(35.0, 139.0) + private val p2 = GeoPoint.fromLatLong(35.0, 140.0) + private val p3 = GeoPoint.fromLatLong(36.0, 140.0) + + @Test + fun circleStrokeWidthDefaultsTo1dp() { + assertEquals(1.dp, CircleState(center = p1, radiusMeters = 100.0).strokeWidth) + } + + @Test + fun polygonStrokeWidthDefaultsTo1dp() { + assertEquals(1.dp, PolygonState(points = listOf(p1, p2, p3)).strokeWidth) + } + + @Test + fun polylineStrokeWidthDefaultsTo1dp() { + assertEquals(1.dp, PolylineState(points = listOf(p1, p2)).strokeWidth) + } + + @Test + fun groundImageOpacityDefaultsTo1() { + val bounds = + GeoRectBounds().apply { + extend(p1) + extend(p3) + } + assertEquals(1.0f, GroundImageState(bounds = bounds, image = FakeDrawable()).opacity) + } + + /** GroundImageState は Drawable を要求するが、テストでは中身を使わない。 */ + private class FakeDrawable : android.graphics.drawable.Drawable() { + override fun draw(canvas: android.graphics.Canvas) = Unit + + override fun setAlpha(alpha: Int) = Unit + + override fun setColorFilter(colorFilter: android.graphics.ColorFilter?) = Unit + + @Deprecated("Deprecated in Java") + override fun getOpacity(): Int = android.graphics.PixelFormat.OPAQUE + } +} diff --git a/src/test/java/com/mapconductor/core/spherical/SignedAreaSignTest.kt b/src/test/java/com/mapconductor/core/spherical/SignedAreaSignTest.kt new file mode 100644 index 0000000..30a1925 --- /dev/null +++ b/src/test/java/com/mapconductor/core/spherical/SignedAreaSignTest.kt @@ -0,0 +1,34 @@ +package com.mapconductor.core.spherical + +import com.mapconductor.core.features.GeoPoint +import com.mapconductor.core.features.GeoPointInterface +import org.junit.Assert.assertEquals +import org.junit.Assert.assertTrue +import org.junit.Test + +/** + * computeSignedArea の符号の約束: 北が上の地図で時計回りのリングが正、反時計回りが負。 + * js-sdk-core の `SignedAreaSign.test.mjs`、ios-sdk の `SignedAreaSignTests` と同じ観点。 + */ +class SignedAreaSignTest { + private val sw = GeoPoint(35.0, 139.0) + private val se = GeoPoint(35.0, 140.0) + private val ne = GeoPoint(36.0, 140.0) + private val nw = GeoPoint(36.0, 139.0) + private val clockwise = listOf(sw, nw, ne, se) + private val counterClockwise = listOf(sw, se, ne, nw) + + private fun check(signedArea: (List) -> Double) { + val cw = signedArea(clockwise) + val ccw = signedArea(counterClockwise) + assertTrue("clockwise should be positive: $cw", cw > 0.0) + assertTrue("counter-clockwise should be negative: $ccw", ccw < 0.0) + assertEquals(cw, -ccw, 1e-3) + } + + @Test + fun spherical_clockwiseIsPositive() = check(Spherical::computeSignedArea) + + @Test + fun wgs84Geodesic_clockwiseIsPositive() = check(WGS84Geodesic::computeSignedArea) +}