[PROJ] About unit conversion
Kristian Evers
kristianevers at gmail.com
Wed May 20 11:28:16 PDT 2020
I’ve look at all the pipeline definitions in the attached code. Below are my comments for each of them.
I agree with Jochem, you seem to be passing linear coordinates to a projection that expects angular
coordinates. Looking back, I see that my comments about the unitconvert steps were not the cause of
your problems. Had you been using different units it could have been the case - it is always better to
be explicit here.
"+proj=utm +lon_0=3w +ellps=GRS80”
No problems here in either directions.
"+proj=pipeline +step +proj=utm +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_out=km
No problems here in either directions, but could benefit from a added +xy_in=m in the unitconvert step.
"+proj=pipeline +step +proj=unitconvert +xy_in=deg +step +proj=utm +lon_0=3w +ellps=GRS80”
Same as above, clarify uniconvert step be adding +xy_out=rad
"+proj=pipeline +step +proj=unitconvert +xy_in=deg +step +proj=utm +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_out=km”
Clarify both unitconvert steps by adding xy_out=rad to the first and +xy_in=m to the second
"+proj=utm +lon_0=3w +ellps=GRS80”
No problems here
"+proj=pipeline +step +proj=unitconvert +xy_in=km +step +proj=utm +lon_0=3w +ellps=GRS80”
Again, clarify unitconvert step by adding xy_out=m
"+proj=pipeline +step +proj=utm +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_out=deg”
This is nonsense, projected coordinates from UTM can’t be converted to degrees.
Probably missing a +inv in the utm step. If that is the case, clarify unitconvert step with +xy_in=rad
"+proj=pipeline +step +proj=unitconvert +xy_in=km +step +proj=utm +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_out=deg”
Also nonsense, passing km to the UTM projection is not possible, UTM want’s radians as input. Probably again missing a +inv
in the utm step. Clarify unitconvert steps by stating both input and output units.
/Kristian
> On 20 May 2020, at 19:30, José Luis García Pallero <jgpallero at gmail.com> wrote:
>
> El mié., 20 may. 2020 a las 19:24, Lesparre, Jochem
> (<Jochem.Lesparre at kadaster.nl <mailto:Jochem.Lesparre at kadaster.nl>>) escribió:
>>
>> I just had a quick look, but to me it looks like you are trying to convert between angles and distances?
>>
>> +proj=utm gives angular results that can be converted to deg, rad or maybe gon, but NOT to m, km of ft! This is such a fundamental part of how PROJ works, that it is probably not explained in the documentation.
>>
>> If I mis understood, sorry for the confusion.
>
> I'm trying to perform first the forward step in the projection, i.e.,
> from geodetic (angular) coordinates to projected (linear) ones. This
> step works without problems whatever the unit combinations I use
> (radians to meters, radians to kilometers, degrees to meters or
> degrees to kilometers). The problem is with the inverse step , i.e.,
> from projected (linear) coordinates to geodetic (angular) ones. I want
> to use the output coordinates of the forward step as input for the
> inverse step. The only combination it works is from meters to radians.
> Whatever other combination fails. Attached to this message
> (https://lists.osgeo.org/pipermail/proj/2020-May/009648.html <https://lists.osgeo.org/pipermail/proj/2020-May/009648.html>) I sent a
> C source code with an example showing what I want to do
>
> Cheers
>
>> Kind regards, Jochem
>>
>>
>> -----Original Message-----
>> From: PROJ <proj-bounces at lists.osgeo.org> On Behalf Of José Luis García Pallero
>> Sent: woensdag 20 mei 2020 18:34
>> To: Kristian Evers <kreve at sdfe.dk>
>> Cc: proj at lists.osgeo.org
>> Subject: Re: [PROJ] About unit conversion
>>
>> El mié., 20 may. 2020 a las 11:00, Kristian Evers (<kreve at sdfe.dk>) escribió:
>>>
>>> I see that you don't specify the input units in your unitconvert
>>> steps. Try adding a
>>> +xy_in parameter to your unitconvert steps and see if that fixes the problem.
>>> When you only specify of of either +xy_in or +xy_out the operation
>>> becomes somewhat uni-directional which is why not all of your pipelines works in reverse.
>>>
>>> See https://proj.org/operations/conversions/unitconvert.html for more.
>>
>> Thank you for your answer, Kristian, but I can't understand why my program fails. Clearly, for FWD and input in radians and output in meters all works without any unitconversion. Also for INV step with meters as input and radians as output, i.e., using only "+proj=utm
>> +lon_0=3w +ellps=GRS80" as parameters string is enough. Also for any
>> combination in the FWD step all works.
>>
>> For example using "+proj=pipeline +step +proj=unitconvert +xy_in=deg
>> +step +proj=utm +lon_0=3w +ellps=GRS80 +step +proj=unitconvert
>> +xy_out=km" accepts as input degrees and returns kilometers. But I
>> cant understand why in this case for the INV step the parameters string "+proj=pipeline +step +proj=unitconvert +xy_in=km +step
>> +proj=utm +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_out=deg"
>> is incorrect. I've tested several combinations and none of them works:
>>
>> "+proj=utm +lon_0=3w +ellps=GRS80" -> wrong results "+proj=pipeline +step +proj=unitconvert +xy_in=km +step +proj=utm
>> +lon_0=3w +ellps=GRS80" -> wrong results
>> "+proj=pipeline +step +proj=unitconvert +xy_in=km +xy_out=m +step
>> +proj=utm +lon_0=3w +ellps=GRS80" -> wrong results (same numeric value
>> as before)
>> "+proj=pipeline +step +proj=unitconvert +xy_in=km +xy_out=deg +step
>> +proj=utm +lon_0=3w +ellps=GRS80" -> proj_create: Error -59
>> (inconsistent unit type between input and output): Pipeline: Bad step
>> definition: proj=unitconvert (inconsistent unit type between input and
>> output)
>> "+proj=pipeline +step +proj=utm +lon_0=3w +ellps=GRS80 +step
>> +proj=unitconvert +xy_in=km +xy_out=deg" -> proj_create: Error -59
>> (inconsistent unit type between input and output): Pipeline: Bad step
>> definition: proj=unitconvert (inconsistent unit type between input and
>> output)
>> "+proj=pipeline +step +proj=unitconvert +xy_out=deg +step +proj=utm
>> +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_in=km" -> wrong
>> results
>> "+proj=pipeline +step +proj=unitconvert +xy_in=km +step +proj=utm
>> +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_out=deg" -> wrong
>> results
>>
>> I can't understand the way as unitconvert works, neither reading https://proj.org/operations/conversions/unitconvert.html
>>
>> For only one conversion (km to radians or meters to degrees) the problem is only again with the INV step
>>
>> Thanks
>>
>>
>>
>>
>>
>>>
>>> /Kristian
>>>
>>>> -----Original Message-----
>>>> From: PROJ <proj-bounces at lists.osgeo.org> On Behalf Of José Luis
>>>> García Pallero
>>>> Sent: 19. maj 2020 10:59
>>>> To: proj at lists.osgeo.org
>>>> Subject: Re: [PROJ] About unit conversion
>>>>
>>>> I've written a small program in order to illustrate my last mail.
>>>> Attached (and pasted at the bottom of the message) I send the code.
>>>> The results are
>>>>
>>>> Original data -> lon: -6.00000, lat: 43.00000 (degrees) -> lon:
>>>> -0.1047197551, lat: 0.7504915784 (radians)
>>>>
>>>> rad to m -> X= 255466.981, Y=4765182.933 (+proj=utm +lon_0=3w
>>>> +ellps=GRS80)
>>>> rad to km -> X= 255.467, Y= 4765.183 (+proj=pipeline +step
>>>> +proj=utm +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_out=km)
>>>> deg to m -> X= 255466.981, Y=4765182.933 (+proj=pipeline +step
>>>> +proj=unitconvert +xy_in=deg +step +proj=utm +lon_0=3w +ellps=GRS80)
>>>> deg to km -> X= 255.467, Y= 4765.183 (+proj=pipeline +step
>>>> +proj=unitconvert +xy_in=deg +step +proj=utm +lon_0=3w +ellps=GRS80
>>>> +step +proj=unitconvert +xy_out=km)
>>>>
>>>> m to rad -> lon=-0.1047197551, lat= 0.7504915784 (+proj=utm
>>>> +lon_0=3w
>>>> +ellps=GRS80)
>>>> km to rad -> lon=-0.0001306633, lat= 0.0000007501 (+proj=pipeline
>>>> +step +proj=unitconvert +xy_in=km +step +proj=utm +lon_0=3w
>>>> +ellps=GRS80)
>>>> m to deg -> lon=-0.1300126568, lat= 0.0130928532 (+proj=pipeline
>>>> +step +proj=utm +lon_0=3w +ellps=GRS80 +step +proj=unitconvert
>>>> +xy_out=deg)
>>>> km to deg -> lon=-0.0001307025, lat= 0.0000000131 (+proj=pipeline
>>>> +step +proj=unitconvert +xy_in=km +step +proj=utm +lon_0=3w
>>>> +ellps=GRS80 +step +proj=unitconvert +xy_out=deg)
>>>>
>>>> As you can see, the only right INV conversion is the
>>>> meters_to_radians one. All other combinations produce all different
>>>> results, but all combinations for the FWD step are right. Im using
>>>> PROJ 7.0.1 from the Debian Sid repositories
>>>>
>>>> ***********************
>>>>
>>>> #include<stdio.h>
>>>> #include<stdlib.h>
>>>> #include<math.h>
>>>> #include<proj.h>
>>>> #define CONST_PI (3.14159265358979323846264338328) int main() {
>>>> double lond1=-6.0,latd1=43.0,lond2=lond1,latd2=latd1;
>>>> double
>>>> lonr1=lond1*CONST_PI/180.0,latr1=latd1*CONST_PI/180.0,lonr2=lonr1,la
>>>> tr2
>>>> =latr1;
>>>> char param_rad_to_m[]="+proj=utm +lon_0=3w +ellps=GRS80";
>>>> char param_rad_to_km[]="+proj=pipeline +step +proj=utm +lon_0=3w
>>>> +ellps=GRS80 +step +proj=unitconvert +xy_out=km";
>>>> char param_deg_to_m[]="+proj=pipeline +step +proj=unitconvert
>>>> +xy_in=deg +step +proj=utm +lon_0=3w +ellps=GRS80";
>>>> char param_deg_to_km[]="+proj=pipeline +step +proj=unitconvert
>>>> +xy_in=deg +step +proj=utm +lon_0=3w +ellps=GRS80 +step
>>>> +proj=unitconvert +xy_out=km";
>>>> char param_m_to_rad[]="+proj=utm +lon_0=3w +ellps=GRS80";
>>>> char param_km_to_rad[]="+proj=pipeline +step +proj=unitconvert
>>>> +xy_in=km +step +proj=utm +lon_0=3w +ellps=GRS80";
>>>> char param_m_to_deg[]="+proj=pipeline +step +proj=utm +lon_0=3w
>>>> +ellps=GRS80 +step +proj=unitconvert +xy_out=deg";
>>>> char param_km_to_deg[]="+proj=pipeline +step +proj=unitconvert
>>>> +xy_in=km +step +proj=utm +lon_0=3w +ellps=GRS80 +step
>>>> +proj=unitconvert +xy_out=deg";
>>>> PJ_CONTEXT *C=proj_context_create();
>>>> PJ *proj=NULL;
>>>> size_t ld=sizeof(double);
>>>>
>>>> //Display original data
>>>> printf("Original data -> lon: %.5lf, lat: %.5lf (degrees) -> lon:
>>>> %.10lf, lat: %.10lf (radians)\n\n",
>>>> lond1,latd1,lonr1,latr1);
>>>> //FWD, radiands to meters
>>>> proj = proj_create(C,param_rad_to_m);
>>>>
>>>> proj_trans_generic(proj,PJ_FWD,&lonr1,ld,1,&latr1,ld,1,NULL,0,0,NULL,0,0);
>>>> printf("rad to m -> X=%11.3lf, Y=%11.3lf
>>>> (%s)\n",lonr1,latr1,param_rad_to_m);
>>>> proj_destroy(proj);
>>>> //FWD, radiands to kilometers
>>>> proj = proj_create(C,param_rad_to_km);
>>>>
>>>> proj_trans_generic(proj,PJ_FWD,&lonr2,ld,1,&latr2,ld,1,NULL,0,0,NULL,0,0);
>>>> printf("rad to km -> X=%11.3lf, Y=%11.3lf
>>>> (%s)\n",lonr2,latr2,param_rad_to_km);
>>>> proj_destroy(proj);
>>>> //FWD, degrees to meters
>>>> proj = proj_create(C,param_deg_to_m);
>>>>
>>>> proj_trans_generic(proj,PJ_FWD,&lond1,ld,1,&latd1,ld,1,NULL,0,0,NULL,0,0);
>>>> printf("deg to m -> X=%11.3lf, Y=%11.3lf
>>>> (%s)\n",lond1,latd1,param_deg_to_m);
>>>> proj_destroy(proj);
>>>> //FWD, degrees to kilometers
>>>> proj = proj_create(C,param_deg_to_km);
>>>>
>>>> proj_trans_generic(proj,PJ_FWD,&lond2,ld,1,&latd2,ld,1,NULL,0,0,NULL,0,0);
>>>> printf("deg to km -> X=%11.3lf, Y=%11.3lf
>>>> (%s)\n\n",lond2,latd2,param_deg_to_km);
>>>> proj_destroy(proj);
>>>> //INV, meters to radians
>>>> proj = proj_create(C,param_m_to_rad);
>>>> proj_trans_generic(proj,PJ_INV,&lonr1,ld,1,&latr1,ld,1,NULL,0,0,NULL,0,0);
>>>> printf("m to rad -> lon=%13.10lf, lat=%13.10lf
>>>> (%s)\n",lonr1,latr1,param_m_to_rad);
>>>> proj_destroy(proj);
>>>> //INV, kilometers to radians
>>>> proj = proj_create(C,param_km_to_rad);
>>>> proj_trans_generic(proj,PJ_INV,&lonr2,ld,1,&latr2,ld,1,NULL,0,0,NULL,0,0);
>>>> printf("km to rad -> lon=%13.10lf, lat=%13.10lf
>>>> (%s)\n",lonr2,latr2,param_km_to_rad);
>>>> proj_destroy(proj);
>>>> //INV, meters to degrees
>>>> proj = proj_create(C,param_m_to_deg);
>>>>
>>>> proj_trans_generic(proj,PJ_INV,&lond1,ld,1,&latd1,ld,1,NULL,0,0,NULL,0,0);
>>>> printf("m to deg -> lon=%13.10lf, lat=%13.10lf
>>>> (%s)\n",lond1,latd1,param_m_to_deg);
>>>> proj_destroy(proj);
>>>> //INV, kilometers to degrees
>>>> proj = proj_create(C,param_km_to_deg);
>>>>
>>>> proj_trans_generic(proj,PJ_INV,&lond2,ld,1,&latd2,ld,1,NULL,0,0,NULL,0,0);
>>>> printf("km to deg -> lon=%13.10lf, lat=%13.10lf
>>>> (%s)\n",lond2,latd2,param_km_to_deg);
>>>> proj_destroy(proj);
>>>> proj_context_destroy(C);
>>>> return 0;
>>>> }
>>>>
>>>> El lun., 18 may. 2020 a las 23:03, José Luis García Pallero
>>>> (<jgpallero at gmail.com>) escribió:
>>>>>
>>>>> Hello:
>>>>>
>>>>> I'm trying to do dome projections (forward and inverse) using the
>>>>> function proj_trans_generic() and parameter string in PROJ format
>>>>> (+
>>>>> style) and I'm a bit confused about unit conversions. My reference
>>>>> for the syntax is
>>>>> https://proj.org/operations/conversions/unitconvert.html
>>>>>
>>>>> For the forward problem I write
>>>>>
>>>>> C = proj_context_create();
>>>>> projection = proj_create(C,param);
>>>>>
>>>> proj_trans_generic(projection,PJ_FWD,lon,sizeof(double),nElem,lat,si
>>>> zeof(d
>>>> ouble),nElem,NULL,0,0,NULL,0,0);
>>>>>
>>>>> I've used the param string
>>>>>
>>>>> '+proj=pipeline +step +proj=unitconvert +xy_in=deg +step +proj=utm
>>>>> +lon_0=3w +ellps=GRS80 +step +proj=unitconvert +xy_out=km'
>>>>>
>>>>> which means that the input coordinates are in degrees and the
>>>>> output in kilometers. All works right. The results are the same if
>>>>> I pass the input coordinates in radians and I use the param string
>>>>>
>>>>> '+proj=pipeline +step +proj=utm +lon_0=3w +ellps=GRS80 +step
>>>>> +proj=unitconvert +xy_out=km'
>>>>>
>>>>> The problem comes when I try to do the inverse step. The procedure
>>>>> is the
>>>> same:
>>>>>
>>>>> C = proj_context_create();
>>>>> projection = proj_create(C,param);
>>>>>
>>>> proj_trans_generic(projection,PJ_INV,x,sizeof(double),nElem,y,sizeof
>>>> (doubl
>>>> e),nElem,NULL,0,0,NULL,0,0);
>>>>>
>>>>> In the input xoordinates are in meters and I use the param string
>>>>>
>>>>> '+proj=utm +lon_0=3w +ellps=GRS80'
>>>>>
>>>>> The resulting longtude and latitude are correct and the units are
>>>>> radians. But if I use the param string
>>>>>
>>>>> '+proj=pipeline +step +proj=utm +lon_0=3w +ellps=GRS80 +step
>>>>> +proj=unitconvert +xy_out=deg'
>>>>>
>>>>> The output is totally wrong. I expect degrees, but the values have
>>>>> no sense. Also, if I use the input X and Y in kilometers ans the
>>>>> param string is
>>>>>
>>>>> '+proj=pipeline +step +proj=unitconvert +xy_in=km +step +proj=utm
>>>>> +lon_0=3w +ellps=GRS80'
>>>>>
>>>>> the result is wrong again.
>>>>>
>>>>> I suppose I'm doing something wrong with unit conversion, but I
>>>>> can not
>>>> find it
>>>>>
>>>>> Can someone help me?
>>>>>
>>>>> Thanks
>>>>>
>>>>> --
>>>>> *****************************************
>>>>> José Luis García Pallero
>>>>> jgpallero at gmail.com
>>>>> (o<
>>>>> / / \
>>>>> V_/_
>>>>> Use Debian GNU/Linux and enjoy!
>>>>> *****************************************
>>>>
>>>>
>>>>
>>>> --
>>>> *****************************************
>>>> José Luis García Pallero
>>>> jgpallero at gmail.com
>>>> (o<
>>>> / / \
>>>> V_/_
>>>> Use Debian GNU/Linux and enjoy!
>>>> *****************************************
>>
>>
>>
>> --
>> *****************************************
>> José Luis García Pallero
>> jgpallero at gmail.com
>> (o<
>> / / \
>> V_/_
>> Use Debian GNU/Linux and enjoy!
>> *****************************************
>> _______________________________________________
>> PROJ mailing list
>> PROJ at lists.osgeo.org
>> https://lists.osgeo.org/mailman/listinfo/proj
>
>
>
> --
> *****************************************
> José Luis García Pallero
> jgpallero at gmail.com <mailto:jgpallero at gmail.com>
> (o<
> / / \
> V_/_
> Use Debian GNU/Linux and enjoy!
> *****************************************
> _______________________________________________
> PROJ mailing list
> PROJ at lists.osgeo.org <mailto:PROJ at lists.osgeo.org>
> https://lists.osgeo.org/mailman/listinfo/proj <https://lists.osgeo.org/mailman/listinfo/proj>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/proj/attachments/20200520/1ea17929/attachment-0001.html>
More information about the PROJ
mailing list