Merge pull request #3503 from rouault/4d_by_ref

Speed-up: pass PJ_COORD by reference for fwd4d and inv4d callbacks

src/4D_api.cpp

@@ -326,8 +326,11 @@ similarly, but prefers the 2D resp. 3D interfaces if available.
                 }
                 P->iCurCoordOp = iBest;
             }
-            PJ_COORD res = direction == PJ_FWD ?
-                        pj_fwd4d( coord, alt.pj ) : pj_inv4d( coord, alt.pj );
+            PJ_COORD res = coord;
+            if( direction == PJ_FWD )
+                pj_fwd4d( res, alt.pj );
+            else
+                pj_inv4d( res, alt.pj );
             if( proj_errno(alt.pj) == PROJ_ERR_OTHER_NETWORK_ERROR ) {
                 return proj_coord_error ();
             }
@@ -368,11 +371,12 @@ similarly, but prefers the 2D resp. 3D interfaces if available.
                         P->iCurCoordOp = i;
                     }
                     if( direction == PJ_FWD ) {
-                        return pj_fwd4d( coord, alt.pj );
+                        pj_fwd4d( coord, alt.pj );
                     }
                     else {
-                        return pj_inv4d( coord, alt.pj );
+                        pj_inv4d( coord, alt.pj );
                     }
+                    return coord;
                 }
             }
         }
@@ -383,9 +387,10 @@ similarly, but prefers the 2D resp. 3D interfaces if available.
 
     P->iCurCoordOp = 0; // dummy value, to be used by proj_trans_get_last_used_operation()
     if (direction == PJ_FWD)
-        return pj_fwd4d (coord, P);
+        pj_fwd4d (coord, P);
     else
-        return pj_inv4d (coord, P);
+        pj_inv4d (coord, P);
+    return coord;
 }
 
 /*****************************************************************************/

src/conversions/axisswap.cpp

@@ -57,6 +57,8 @@ It is only necessary to specify the axes that are affected by the swap
 #include <stdlib.h>
 #include <string.h>
 
+#include <algorithm>
+
 #include "proj.h"
 #include "proj_internal.h"
 
@@ -75,20 +77,14 @@ static int sign(int x) {
 
 static PJ_XY forward_2d(PJ_LP lp, PJ *P) {
     struct pj_opaque *Q = (struct pj_opaque *) P->opaque;
-    unsigned int i;
-    PJ_COORD out, in;
-
-    in.v[0] = lp.lam;
-    in.v[1] = lp.phi;
-    out = proj_coord_error();
+    PJ_XY xy;
 
-    for (i=0; i<2; i++)
-        out.v[i] = in.v[Q->axis[i]] * Q->sign[i];
-
-    return out.xy;
+    double in[2] = {lp.lam, lp.phi};
+    xy.x = in[Q->axis[0]] * Q->sign[0];
+    xy.y = in[Q->axis[1]] * Q->sign[1];
+    return xy;
 }
 
-
 static PJ_LP reverse_2d(PJ_XY xy, PJ *P) {
     struct pj_opaque *Q = (struct pj_opaque *) P->opaque;
     unsigned int i;
@@ -104,7 +100,6 @@ static PJ_LP reverse_2d(PJ_XY xy, PJ *P) {
     return out.lp;
 }
 
-
 static PJ_XYZ forward_3d(PJ_LPZ lpz, PJ *P) {
     struct pj_opaque *Q = (struct pj_opaque *) P->opaque;
     unsigned int i;
@@ -137,32 +132,31 @@ static PJ_LPZ reverse_3d(PJ_XYZ xyz, PJ *P) {
     return out.lpz;
 }
 
+static void swap_xy_4d(PJ_COORD& coo, PJ *) {
+    std::swap(coo.xyzt.x, coo.xyzt.y);
+}
+
 
-static PJ_COORD forward_4d(PJ_COORD coo, PJ *P) {
+static void forward_4d(PJ_COORD& coo, PJ *P) {
     struct pj_opaque *Q = (struct pj_opaque *) P->opaque;
     unsigned int i;
     PJ_COORD out;
 
-    out = proj_coord_error();
-
     for (i=0; i<4; i++)
         out.v[i] = coo.v[Q->axis[i]] * Q->sign[i];
-
-    return out;
+    coo = out;
 }
 
 
-static PJ_COORD reverse_4d(PJ_COORD coo, PJ *P) {
+static void reverse_4d(PJ_COORD& coo, PJ *P) {
     struct pj_opaque *Q = (struct pj_opaque *) P->opaque;
     unsigned int i;
     PJ_COORD out;
 
-    out = proj_coord_error();
-
     for (i=0; i<4; i++)
         out.v[Q->axis[i]] = coo.v[i] * Q->sign[i];
 
-    return out;
+    coo = out;
 }
 
 
@@ -278,9 +272,18 @@ PJ *CONVERSION(axisswap,0) {
         P->fwd3d  = forward_3d;
         P->inv3d  = reverse_3d;
     }
-    if (n == 2 && Q->axis[0] < 2 && Q->axis[1] < 2) {
-        P->fwd    = forward_2d;
-        P->inv    = reverse_2d;
+    if (n == 2) {
+        if( Q->axis[0] == 1 && Q->sign[0] == 1 &&
+            Q->axis[1] == 0 && Q->sign[1] == 1 )
+        {
+            P->fwd4d = swap_xy_4d;
+            P->inv4d = swap_xy_4d;
+        }
+        else if( Q->axis[0] < 2 && Q->axis[1] < 2 )
+        {
+            P->fwd    = forward_2d;
+            P->inv    = reverse_2d;
+        }
     }
 
 

src/conversions/geoc.cpp

@@ -36,13 +36,13 @@
 PROJ_HEAD(geoc, "Geocentric Latitude");
 
 /* Geographical to geocentric */
-static PJ_COORD forward(PJ_COORD coo, PJ *P) {
-    return pj_geocentric_latitude (P, PJ_FWD, coo);
+static void forward(PJ_COORD& coo, PJ *P) {
+    coo = pj_geocentric_latitude (P, PJ_FWD, coo);
 }
 
 /* Geocentric to geographical */
-static PJ_COORD inverse(PJ_COORD coo, PJ *P) {
-    return pj_geocentric_latitude (P, PJ_INV, coo);
+static void inverse(PJ_COORD& coo, PJ *P) {
+    coo = pj_geocentric_latitude (P, PJ_INV, coo);
 }
 
 

src/conversions/noop.cpp

@@ -4,9 +4,7 @@
 
 PROJ_HEAD(noop,    "No operation");
 
-static PJ_COORD noop(PJ_COORD coord, PJ *P) {
-    (void) P;
-    return coord;
+static void noop(PJ_COORD& , PJ *) {
 }
 
 PJ *CONVERSION(noop, 0) {

src/conversions/set.cpp

@@ -19,7 +19,7 @@ struct Set {
 };
 } // anonymous namespace
 
-static PJ_COORD set_fwd_inv(PJ_COORD point, PJ *P) {
+static void set_fwd_inv(PJ_COORD& point, PJ *P) {
 
     struct Set *set = static_cast<struct Set*>(P->opaque);
 
@@ -31,8 +31,6 @@ static PJ_COORD set_fwd_inv(PJ_COORD point, PJ *P) {
         point.v[2] = set->v3_val;
     if (set->v4)
         point.v[3] = set->v4_val;
-
-    return point;
 }
 
 PJ *OPERATION(set, 0) {

src/conversions/topocentric.cpp

@@ -48,28 +48,27 @@ struct pj_opaque {
 } // anonymous namespace
 
 // Convert from geocentric to topocentric
-static PJ_COORD topocentric_fwd(PJ_COORD in, PJ * P)
+static void topocentric_fwd(PJ_COORD& coo, PJ * P)
 {
     struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque);
-    PJ_COORD out;
-    const double dX = in.xyz.x - Q->X0;
-    const double dY = in.xyz.y - Q->Y0;
-    const double dZ = in.xyz.z - Q->Z0;
-    out.xyz.x = -dX * Q->sinlam0              + dY * Q->coslam0;
-    out.xyz.y = -dX * Q->sinphi0 * Q->coslam0 - dY * Q->sinphi0 * Q->sinlam0 + dZ * Q->cosphi0;
-    out.xyz.z =  dX * Q->cosphi0 * Q->coslam0 + dY * Q->cosphi0 * Q->sinlam0 + dZ * Q->sinphi0;
-    return out;
+    const double dX = coo.xyz.x - Q->X0;
+    const double dY = coo.xyz.y - Q->Y0;
+    const double dZ = coo.xyz.z - Q->Z0;
+    coo.xyz.x = -dX * Q->sinlam0              + dY * Q->coslam0;
+    coo.xyz.y = -dX * Q->sinphi0 * Q->coslam0 - dY * Q->sinphi0 * Q->sinlam0 + dZ * Q->cosphi0;
+    coo.xyz.z =  dX * Q->cosphi0 * Q->coslam0 + dY * Q->cosphi0 * Q->sinlam0 + dZ * Q->sinphi0;
 }
 
 // Convert from topocentric to geocentric
-static PJ_COORD topocentric_inv(PJ_COORD in, PJ * P)
+static void topocentric_inv(PJ_COORD& coo, PJ * P)
 {
     struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque);
-    PJ_COORD out;
-    out.xyz.x = Q->X0 - in.xyz.x * Q->sinlam0 - in.xyz.y * Q->sinphi0 * Q->coslam0 + in.xyz.z * Q->cosphi0 * Q->coslam0;
-    out.xyz.y = Q->Y0 + in.xyz.x * Q->coslam0 - in.xyz.y * Q->sinphi0 * Q->sinlam0 + in.xyz.z * Q->cosphi0 * Q->sinlam0;
-    out.xyz.z = Q->Z0                         + in.xyz.y * Q->cosphi0              + in.xyz.z * Q->sinphi0;
-    return out;
+    const double x = coo.xyz.x;
+    const double y = coo.xyz.y;
+    const double z = coo.xyz.z;
+    coo.xyz.x = Q->X0 - x * Q->sinlam0 - y * Q->sinphi0 * Q->coslam0 + z * Q->cosphi0 * Q->coslam0;
+    coo.xyz.y = Q->Y0 + x * Q->coslam0 - y * Q->sinphi0 * Q->sinlam0 + z * Q->cosphi0 * Q->sinlam0;
+    coo.xyz.z = Q->Z0                         + y * Q->cosphi0              + z * Q->sinphi0;
 }
 
 

src/conversions/unitconvert.cpp

@@ -351,42 +351,36 @@ static PJ_LPZ reverse_3d(PJ_XYZ xyz, PJ *P) {
 
 
 /***********************************************************************/
-static PJ_COORD forward_4d(PJ_COORD obs, PJ *P) {
+static void forward_4d(PJ_COORD& coo, PJ *P) {
 /************************************************************************
     Forward conversion of time units
 ************************************************************************/
     struct pj_opaque_unitconvert *Q = (struct pj_opaque_unitconvert *) P->opaque;
-    PJ_COORD out = obs;
 
     /* delegate unit conversion of physical dimensions to the 3D function */
-    out.xyz = forward_3d(obs.lpz, P);
+    coo.xyz = forward_3d(coo.lpz, P);
 
     if (Q->t_in_id >= 0)
-        out.xyzt.t = time_units[Q->t_in_id].t_in( obs.xyzt.t );
+        coo.xyzt.t = time_units[Q->t_in_id].t_in( coo.xyzt.t );
     if (Q->t_out_id >= 0)
-        out.xyzt.t = time_units[Q->t_out_id].t_out( out.xyzt.t );
-
-    return out;
+        coo.xyzt.t = time_units[Q->t_out_id].t_out( coo.xyzt.t );
 }
 
 
 /***********************************************************************/
-static PJ_COORD reverse_4d(PJ_COORD obs, PJ *P) {
+static void reverse_4d(PJ_COORD& coo, PJ *P) {
 /************************************************************************
     Reverse conversion of time units
 ************************************************************************/
     struct pj_opaque_unitconvert *Q = (struct pj_opaque_unitconvert *) P->opaque;
-    PJ_COORD out = obs;
 
     /* delegate unit conversion of physical dimensions to the 3D function */
-    out.lpz = reverse_3d(obs.xyz, P);
+    coo.lpz = reverse_3d(coo.xyz, P);
 
     if (Q->t_out_id >= 0)
-        out.xyzt.t = time_units[Q->t_out_id].t_in( obs.xyzt.t );
+        coo.xyzt.t = time_units[Q->t_out_id].t_in( coo.xyzt.t );
     if (Q->t_in_id >= 0)
-        out.xyzt.t = time_units[Q->t_in_id].t_out( out.xyzt.t );
-
-    return out;
+        coo.xyzt.t = time_units[Q->t_in_id].t_out( coo.xyzt.t );
 }
 
 /***********************************************************************/

src/fwd.cpp

@@ -203,7 +203,7 @@ PJ_XY pj_fwd(PJ_LP lp, PJ *P) {
         coo.xyz = xyz;
     }
     else if (P->fwd4d)
-        coo = P->fwd4d (coo, P);
+        P->fwd4d (coo, P);
     else {
         proj_errno_set (P, PROJ_ERR_OTHER_NO_INVERSE_OP);
         return proj_coord_error ().xy;
@@ -238,7 +238,7 @@ PJ_XYZ pj_fwd3d(PJ_LPZ lpz, PJ *P) {
         coo.xyz = xyz;
     }
     else if (P->fwd4d)
-        coo = P->fwd4d (coo, P);
+        P->fwd4d (coo, P);
     else if (P->fwd)
     {
         const auto xy = P->fwd (coo.lp, P);
@@ -259,19 +259,22 @@ PJ_XYZ pj_fwd3d(PJ_LPZ lpz, PJ *P) {
 
 
 
-PJ_COORD pj_fwd4d (PJ_COORD coo, PJ *P) {
+bool pj_fwd4d (PJ_COORD& coo, PJ *P) {
 
     const int last_errno = P->ctx->last_errno;
     P->ctx->last_errno = 0;
 
     if (!P->skip_fwd_prepare)
         fwd_prepare (P, coo);
     if (HUGE_VAL==coo.v[0])
-        return proj_coord_error ();
+    {
+        coo = proj_coord_error ();
+        return false;
+    }
 
     /* Call the highest dimensional converter available */
     if (P->fwd4d)
-        coo = P->fwd4d (coo, P);
+        P->fwd4d (coo, P);
     else if (P->fwd3d)
     {
         const auto xyz = P->fwd3d (coo.lpz, P);
@@ -284,13 +287,24 @@ PJ_COORD pj_fwd4d (PJ_COORD coo, PJ *P) {
     }
     else {
         proj_errno_set (P, PROJ_ERR_OTHER_NO_INVERSE_OP);
-        return proj_coord_error ();
+        coo = proj_coord_error ();
+        return false;
     }
     if (HUGE_VAL==coo.v[0])
-        return proj_coord_error ();
+    {
+        coo = proj_coord_error ();
+        return false;
+    }
 
     if (!P->skip_fwd_finalize)
         fwd_finalize (P, coo);
 
-    return error_or_coord(P, coo, last_errno);
+    if (P->ctx->last_errno)
+    {
+        coo = proj_coord_error();
+        return false;
+    }
+
+    P->ctx->last_errno = last_errno;
+    return true;
 }