[PROJ] Angular units in PROJ for +proj=latlong when using proj_create_crs_to_crs() and proj_create()

José Luis García Pallero jgpallero at gmail.com
Sun Feb 23 10:45:15 PST 2020


Hello:

Attached I send a program that uses pj_trans() in order to perform a
crs to crs transformation and a cartographic projection which are both
the same operation. The projection is defined as "+proj=utm +lon_0=3w
+ellps=GRS80", and the crs to crs as origin: "+proj=longlat
+ellps=GRS80", destination "+proj=utm +lon_0=3w +ellps=GRS80". In both
cases the PJ_FWD step is used

But I found that when I use pj_trans() as projection the geodetic
coordinates must be given in radians, but when the crs2crs operation
is performed the angular units must be degrees. The PROJ documentation
says (https://proj.org/usage/projections.html?highlight=latlong):

"
A few special projections deviate from this behaviour, most notably
the latlong pseudo-projection that returns degrees.
"

But as it can be seen in the attached code this behavior is only in
the crs2crs operation, so I suppose the documentation at this point is
not exact. But also the "most notably" text suggest that there are
other definitions that works in degrees? What is the reason for this
"inconsistency"? Wouldn't it be more practical to work always in
radians, at least at API level?

Best regards

#include<stdio.h>
#include<math.h>
#include<proj.h>
main(int argc, char **argv)
{
    PJ *P1,*P2;
    PJ_COORD c1,c2,c1_out,c2_out;
    P1 = proj_create_crs_to_crs(PJ_DEFAULT_CTX,
                               "+proj=longlat +ellps=GRS80",
                               "+proj=utm +lon_0=3w +ellps=GRS80",
                               NULL);
    P2 = proj_create(PJ_DEFAULT_CTX,"+proj=utm +lon_0=3w +ellps=GRS80");
    PJ* P_for_GIS = proj_normalize_for_visualization(PJ_DEFAULT_CTX,P1);
    proj_destroy(P1);
    P1 = P_for_GIS;
    c1.lpz.lam = -3.0;
    c1.lpz.phi = 45.0;
    c1.lpz.z = 100.0;
    c2.lpz.lam = -3.0*M_PI/180.0;
    c2.lpz.phi = 45.0*M_PI/180.0;
    c2.lpz.z = 100.0;
    c1_out = proj_trans(P1,PJ_FWD,c1);
    c2_out = proj_trans(P2,PJ_FWD,c2);
    printf("P1: x=%.3f y=%.3f z=%.3f\n",
           c1_out.xyz.x,c1_out.xyz.y,c1_out.xyz.z);
    printf("P2: x=%.3f y=%.3f z=%.3f\n",
           c2_out.xyz.x,c2_out.xyz.y,c2_out.xyz.z);
    proj_destroy(P1);
    proj_destroy(P2);
    return 0;
}

-- 
*****************************************
José Luis García Pallero
jgpallero at gmail.com
(o<
/ / \
V_/_
Use Debian GNU/Linux and enjoy!
*****************************************
-------------- next part --------------
A non-text attachment was scrubbed...
Name: test_proj_units.c
Type: text/x-csrc
Size: 1016 bytes
Desc: not available
URL: <http://lists.osgeo.org/pipermail/proj/attachments/20200223/9e2cfb16/attachment-0001.c>


More information about the PROJ mailing list