<div dir="ltr"><div>About proj_create_crs_to_crs giving wrong results, I do not know. Maybe an axis order problem?</div><div>Providing a transformation instead of a CRS is obviously faster. Less work for the algorithm to obtain the transformation from the CRS (creating temporary structures, etc). I just mentioned the option to provide a CRS for completeness, not expecting that it is the fastest solution.</div><div>As Thomas mentioned, internally the algorithm is using a projection transformation.</div><div><br></div><div>Cheers</div><div>Javier<br></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Fri, 25 Oct 2024 at 08:28, Roger Oberholtzer <<a href="mailto:roger.oberholtzer@gmail.com">roger.oberholtzer@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">On Thu, Oct 24, 2024 at 5:42 PM Javier Jimenez Shaw <<a href="mailto:j1@jimenezshaw.com" target="_blank">j1@jimenezshaw.com</a>> wrote:<br>
><br>
> <a href="https://proj.org/en/9.5/development/reference/functions.html#c.proj_factors" rel="noreferrer" target="_blank">https://proj.org/en/9.5/development/reference/functions.html#c.proj_factors</a><br>
> You just provide the projected CRS, just created from the EPSG code. It converts it internally into the projection over its geographic.<br>
<br>
If I use something defined with proj_create_crs_to_crs I get incorrect<br>
values for the meridian convergence. If I provide an EPSG code to<br>
proj_create and use that, I get the correct values, but it takes an<br>
extremely long time for each calculation. As I wrote, projections that<br>
should take 15 minutes take 15 hours. The thing that seems to work<br>
best is proj_create with a proj string that does not contain<br>
+type=crs. I am wondering if this is all expected behavior.<br>
<br>
<br>
><br>
> On Thu, 24 Oct 2024 at 15:19, Roger Oberholtzer <<a href="mailto:roger.oberholtzer@gmail.com" target="_blank">roger.oberholtzer@gmail.com</a>> wrote:<br>
>><br>
>> On Thu, Oct 24, 2024 at 2:36 PM Javier Jimenez Shaw <<a href="mailto:j1@jimenezshaw.com" target="_blank">j1@jimenezshaw.com</a>> wrote:<br>
>> ><br>
>> > The last versions also accept a projected CRS, extracting the projection from it.<br>
>><br>
>> I am using 9.4.0. Do you mean that I should be able to use this?<br>
>><br>
>> proj_normalize_for_visualization(0,<br>
>>                                        proj_create_crs_to_crs( 0,<br>
>>                                                                "EPSG:4326",<br>
>>                                                                "EPSG:3006",<br>
>>                                                                0))))<br>
>><br>
>> If that's what you mean, it seems that the convergence values are not<br>
>> as expected. Or do you mean something else?<br>
>><br>
>> ><br>
>> > On Thu, 24 Oct 2024, 14:27 Thomas Knudsen via PROJ, <<a href="mailto:proj@lists.osgeo.org" target="_blank">proj@lists.osgeo.org</a>> wrote:<br>
>> >><br>
>> >> Roger,<br>
>> >><br>
>> >> The "factors" elements refer to a projection, not to a CRS<br>
>> >> (the CRS includes reference frame information), and not<br>
>> >> to a transformation (which describes a path between two<br>
>> >> different CRS).<br>
>> >><br>
>> >> So you should only call proj_factors on a PJ object obtained<br>
>> >> from a plain classic proj(ection) string, i.e:<br>
>> >><br>
>> >> p = proj_create(0, "proj=utm  zone=33  ellps=GRS80");<br>
>> >><br>
>> >><br>
>> >> Den tors. 24. okt. 2024 kl. 14.11 skrev Roger Oberholtzer via PROJ <<a href="mailto:proj@lists.osgeo.org" target="_blank">proj@lists.osgeo.org</a>>:<br>
>> >>><br>
>> >>> Hi all,<br>
>> >>><br>
>> >>> I have been using proj_factors in a couple CRS to get the<br>
>> >>> median_convergence value. I need this when rotating a point that has<br>
>> >>> been projected into the CRS.<br>
>> >>><br>
>> >>> In my usual setup, I define something like the following (WGS84 (from<br>
>> >>> a GPS receiver) to Sweref99):<br>
>> >>><br>
>> >>> proj_normalize_for_visualization(0,<br>
>> >>>                                        proj_create_crs_to_crs( 0,<br>
>> >>>                                                                "EPSG:4326",<br>
>> >>>                                                                "EPSG:3006",<br>
>> >>>                                                                0))))<br>
>> >>> return(ISERROR);<br>
>> >>><br>
>> >>> This works fine. I get the points I expect. But as I will need to<br>
>> >>> rotate that point using data from the GPS (heading is true north), I<br>
>> >>> need to know how that differs from the projected point's north. So I<br>
>> >>> use proj_factors. As it turns out, I cannot use the definition from<br>
>> >>> above to do this. The median_convergence value is wrong. So what I<br>
>> >>> have been doing (for EPSG:3006 - it is CRS-specific) is defining<br>
>> >>> something like this:<br>
>> >>><br>
>> >>> proj_create(0, "+proj=utm +zone=33 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0<br>
>> >>> +units=m +no_defs");<br>
>> >>><br>
>> >>> You will note that this is the implementation of EPSG:3006, but<br>
>> >>> without +type=crs. Why without that parameter? Time. When including<br>
>> >>> it, the proj_factors call takes around 60 times longer. This is<br>
>> >>> extremely significant when generating, say, point clouds with hundreds<br>
>> >>> of million points. A processing time of 15 minutes becomes 15 hours.<br>
>> >>> The values provided are exactly the same. But one takes significantly<br>
>> >>> longer to provide them.<br>
>> >>><br>
>> >>> I see this effect in more than EPSG:3006. For example EPSG:5677. I<br>
>> >>> suspect that it is generally the case.<br>
>> >>><br>
>> >>> What I would like to do is something like this:<br>
>> >>><br>
>> >>> proj_normalize_for_visualization(0,<br>
>> >>>                                        proj_create_crs_to_crs( 0,<br>
>> >>>                                                                "EPSG:4326",<br>
>> >>>                                                                "EPSG:3006",<br>
>> >>>                                                                0))))<br>
>> >>> return(ISERROR);<br>
>> >>><br>
>> >>> proj_create(0, "EPSG:3006");<br>
>> >>><br>
>> >>><br>
>> >>> Where the first is used to project the point, and the second is used<br>
>> >>> to call proj_factors for that point.<br>
>> >>><br>
>> >>> This would allow a general solution, which is what I am after. So I<br>
>> >>> either need to call something other than proj_create, or I need to see<br>
>> >>> how to remove +type=crs in this context.<br>
>> >>><br>
>> >>> Does this make sense?<br>
>> >>><br>
>> >>> --<br>
>> >>> Roger Oberholtzer<br>
>> >>> _______________________________________________<br>
>> >>> PROJ mailing list<br>
>> >>> <a href="mailto:PROJ@lists.osgeo.org" target="_blank">PROJ@lists.osgeo.org</a><br>
>> >>> <a href="https://lists.osgeo.org/mailman/listinfo/proj" rel="noreferrer" target="_blank">https://lists.osgeo.org/mailman/listinfo/proj</a><br>
>> >><br>
>> >> _______________________________________________<br>
>> >> PROJ mailing list<br>
>> >> <a href="mailto:PROJ@lists.osgeo.org" target="_blank">PROJ@lists.osgeo.org</a><br>
>> >> <a href="https://lists.osgeo.org/mailman/listinfo/proj" rel="noreferrer" target="_blank">https://lists.osgeo.org/mailman/listinfo/proj</a><br>
>><br>
>><br>
>><br>
>> --<br>
>> Roger Oberholtzer<br>
<br>
<br>
<br>
-- <br>
Roger Oberholtzer<br>
</blockquote></div></div>