|
| 1 | +vec3 projectToGeographic(vec3 cartographic) |
| 2 | +{ |
| 3 | + float lon = cartographic.y; |
| 4 | + float lat = cartographic.x; |
| 5 | + float altitude = cartographic.z; |
| 6 | + |
| 7 | + float plate_x = u_ellipsoidRadii.x * lat; |
| 8 | + float plate_y = u_ellipsoidRadii.x * lon; |
| 9 | + |
| 10 | + return vec3(altitude, plate_x, plate_y); |
| 11 | +} |
| 12 | + |
| 13 | +vec3 cartesianToCartographic(vec3 cartesian) { |
| 14 | + // Correct the incoming (y, z, x) coordinate order to standard ECEF (x, y, z). |
| 15 | + vec3 ecef = vec3(cartesian.z, cartesian.x, cartesian.y); |
| 16 | + |
| 17 | + // Check if the ellipsoid is a sphere, which allows for a much faster, direct calculation. |
| 18 | + bool isSphere = (u_ellipsoidRadii.x == u_ellipsoidRadii.y) && (u_ellipsoidRadii.x == u_ellipsoidRadii.z); |
| 19 | + |
| 20 | + if (isSphere) { |
| 21 | + float mag = length(ecef); |
| 22 | + // Avoid division by zero at the center |
| 23 | + if (mag < czm_epsilon6) { |
| 24 | + return vec3(0.0, 0.0, -u_ellipsoidRadii.x); |
| 25 | + } |
| 26 | + float longitude = atan(ecef.y, ecef.x); |
| 27 | + float latitude = asin(ecef.z / mag); |
| 28 | + float height = mag - u_ellipsoidRadii.x; |
| 29 | + return vec3(longitude, latitude, height); |
| 30 | + } |
| 31 | + |
| 32 | + // --- Ellipsoidal Calculation (Iterative Method) --- |
| 33 | + |
| 34 | + // The vector from the center to the point on the XY plane |
| 35 | + vec2 p = ecef.xy; |
| 36 | + float pLength = length(p); |
| 37 | + |
| 38 | + // If the point is near the Z-axis, the calculation is straightforward |
| 39 | + if (pLength < czm_epsilon6) { |
| 40 | + // Longitude is undefined at the poles, but 0.0 is a standard convention. |
| 41 | + float longitude = 0.0; |
| 42 | + float latitude = (ecef.z >= 0.0) ? czm_pi / 2.0 : -czm_pi / 2.0; |
| 43 | + float height = abs(ecef.z) - u_ellipsoidRadii.z; |
| 44 | + return vec3(longitude, latitude, height); |
| 45 | + } |
| 46 | + |
| 47 | + // Key ellipsoid parameters |
| 48 | + float a = u_ellipsoidRadii.x; |
| 49 | + float c = u_ellipsoidRadii.z; |
| 50 | + float a2 = u_ellipsoidRadiiSquared.x; |
| 51 | + float c2 = u_ellipsoidRadiiSquared.z; |
| 52 | + float e2 = (a2 - c2) / a2; // First eccentricity squared |
| 53 | + |
| 54 | + // Initial guess for latitude is the geocentric latitude |
| 55 | + float latitude = atan(ecef.z, pLength); |
| 56 | + float N; // Radius of curvature in the prime vertical |
| 57 | + |
| 58 | + // Iteratively refine the latitude. 3-4 iterations are sufficient for high precision. |
| 59 | + for (int i = 0; i < 4; ++i) { |
| 60 | + float sinLat = sin(latitude); |
| 61 | + N = a / sqrt(1.0 - e2 * sinLat * sinLat); |
| 62 | + latitude = atan((ecef.z + N * e2 * sinLat) / pLength); |
| 63 | + } |
| 64 | + |
| 65 | + // With the final latitude, calculate longitude and height |
| 66 | + float longitude = atan(ecef.y, ecef.x); |
| 67 | + float sinLat = sin(latitude); |
| 68 | + float cosLat = cos(latitude); |
| 69 | + N = a / sqrt(1.0 - e2 * sinLat * sinLat); |
| 70 | + float height = (pLength / cosLat) - N; |
| 71 | + |
| 72 | + return vec3(longitude, latitude, height); |
| 73 | +} |
| 74 | + |
1 | 75 | vec4 geometryStage(inout ProcessedAttributes attributes, mat4 modelView, mat3 normal)
|
2 | 76 | {
|
3 | 77 | vec4 computedPosition;
|
4 | 78 |
|
5 | 79 | // Compute positions in different coordinate systems
|
6 | 80 | vec3 positionMC = attributes.positionMC;
|
7 | 81 | v_positionMC = positionMC;
|
8 |
| - v_positionEC = (modelView * vec4(positionMC, 1.0)).xyz; |
9 | 82 |
|
10 | 83 | #if defined(USE_2D_POSITIONS) || defined(USE_2D_INSTANCING)
|
11 | 84 | vec3 position2D = attributes.position2D;
|
12 | 85 | vec3 positionEC = (u_modelView2D * vec4(position2D, 1.0)).xyz;
|
13 | 86 | computedPosition = czm_projection * vec4(positionEC, 1.0);
|
14 | 87 | #else
|
15 |
| - computedPosition = czm_projection * vec4(v_positionEC, 1.0); |
| 88 | + if(czm_morphTime != 1.){ |
| 89 | + vec3 cartographic = cartesianToCartographic(positionMC); |
| 90 | + vec3 position = projectToGeographic(cartographic); |
| 91 | + |
| 92 | + vec3 positionEC = (czm_view * vec4(position, 1.0)).xyz; |
| 93 | + computedPosition = czm_projection * vec4(positionEC, 1.0); |
| 94 | + }else{ |
| 95 | + v_positionEC = (modelView * vec4(positionMC, 1.0)).xyz; |
| 96 | + computedPosition = czm_projection * vec4(v_positionEC, 1.0); |
| 97 | + } |
16 | 98 | #endif
|
17 | 99 |
|
18 | 100 | // Sometimes the custom shader and/or style needs this
|
|
0 commit comments