<html xmlns:v="urn:schemas-microsoft-com:vml" xmlns:o="urn:schemas-microsoft-com:office:office" xmlns:w="urn:schemas-microsoft-com:office:word" xmlns:m="http://schemas.microsoft.com/office/2004/12/omml" xmlns="http://www.w3.org/TR/REC-html40">
<head>
<meta http-equiv="Content-Type" content="text/html; charset=utf-8">
<meta name="Generator" content="Microsoft Word 15 (filtered medium)">
<style><!--
/* Font Definitions */
@font-face
{font-family:"Cambria Math";
panose-1:2 4 5 3 5 4 6 3 2 4;}
@font-face
{font-family:Calibri;
panose-1:2 15 5 2 2 2 4 3 2 4;}
@font-face
{font-family:Consolas;
panose-1:2 11 6 9 2 2 4 3 2 4;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
{margin:0in;
font-size:11.0pt;
font-family:"Calibri",sans-serif;
mso-ligatures:standardcontextual;}
a:link, span.MsoHyperlink
{mso-style-priority:99;
color:#0563C1;
text-decoration:underline;}
pre
{mso-style-priority:99;
mso-style-link:"HTML Preformatted Char";
margin:0in;
font-size:10.0pt;
font-family:"Courier New";
mso-ligatures:none;}
span.HTMLPreformattedChar
{mso-style-name:"HTML Preformatted Char";
mso-style-priority:99;
mso-style-link:"HTML Preformatted";
font-family:"Consolas",serif;
mso-ligatures:standardcontextual;}
span.EmailStyle22
{mso-style-type:personal-reply;
font-family:"Calibri",sans-serif;
color:windowtext;}
span.pre
{mso-style-name:pre;}
.MsoChpDefault
{mso-style-type:export-only;
font-size:10.0pt;}
@page WordSection1
{size:8.5in 11.0in;
margin:1.0in 1.0in 1.0in 1.0in;}
div.WordSection1
{page:WordSection1;}
--></style><!--[if gte mso 9]><xml>
<o:shapedefaults v:ext="edit" spidmax="1026" />
</xml><![endif]--><!--[if gte mso 9]><xml>
<o:shapelayout v:ext="edit">
<o:idmap v:ext="edit" data="1" />
</o:shapelayout></xml><![endif]-->
</head>
<body lang="EN-US" link="#0563C1" vlink="purple" style="word-wrap:break-word">
<div class="WordSection1">
<p class="MsoNormal">Eve, thanks for looking to this but I think I’ve not been clear enough in my goals.<o:p></o:p></p>
<p class="MsoNormal">I don’t want to clip into the [0 360] interval. Quite on the contrary, I want to make maps that extend that limit but keep ofc the earth periodicity. My initial attempts worked quite well with the Van Der Grinten projection (<a href="https://forum.generic-mapping-tools.org/t/best-projection-for-rectangular-world-map/3715/82?u=joaquim">https://forum.generic-mapping-tools.org/t/best-projection-for-rectangular-world-map/3715/82?u=joaquim</a>).<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">For the Winkel Triple gdalwarp actually works quite well too<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">try this with the grid in the issue <a href="https://github.com/OSGeo/gdal/issues/7644">
https://github.com/OSGeo/gdal/issues/7644</a><o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">gdalwarp -t_srs "+proj=wintri +over" earth15.grd lixo.tiff<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">It’s the vector step that fails (but the cl540 has 384 segments, not just one).<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">Anyway, this is not an important matter. It would just be nice if I could to with “+proj=wintry +over” what I can with “+proj=vandg +over”<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">Note, I saw the warning at <a href="https://proj.org/usage/projections.html#longitude-wrapping">
https://proj.org/usage/projections.html#longitude-wrapping</a><o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">“Note however that for most projections where the 180 meridian does not project to a straight line,
<span class="pre"><span style="font-size:10.0pt;font-family:"Courier New"">+over</span></span> will have no effect or will not lead to expected results.”<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">The Van Der Grinten doesn’t fulfill that condition but things still work for it, so I don’t really know how to interpret that warning.<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">Joaquim<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<div>
<div style="border:none;border-top:solid #E1E1E1 1.0pt;padding:3.0pt 0in 0in 0in">
<p class="MsoNormal"><b><span style="mso-ligatures:none">From:</span></b><span style="mso-ligatures:none"> Even Rouault <even.rouault@spatialys.com>
<br>
<b>Sent:</b> Sunday, April 30, 2023 2:43 PM<br>
<b>To:</b> Joaquim Manuel Freire Luís <jluis@ualg.pt>; gdal-dev@lists.osgeo.org<br>
<b>Subject:</b> Re: [gdal-dev] +proj=wintri +over option works weel in gdalwarp but not in ogr2ogr<o:p></o:p></span></p>
</div>
</div>
<p class="MsoNormal"><o:p> </o:p></p>
<p>Joaquim,<span style="mso-ligatures:none"><o:p></o:p></span></p>
<p>I've looked at that. Experimentally, it doesn't seem that the current maths for wintri in +over mode give reasonable results for longitudes > 360deg in absolute value (and I'm not sure if there would be a natural extension possible. Far beyond my knowledge
of the underlying maths).<o:p></o:p></p>
<p>gdalwarp does a little better by default because it couples the forward and inverse transformer to check the domain of validity and experimentally must find this limitation and will thus give up for source pixels whose longitude > 360 deg. For ogr2ogr /
vector this is more complicated as reprojection of a single geometry (here cl540.gpkg is a single polygon with longitudes in [-540,540] range). is a "all points succeed to reproject => OK" / "one single point fails to reproject => KO" type of situation, at
least in the logic followed by OGRGeometry::transform() logic.<o:p></o:p></p>
<p>So while the following patch to PROJ is probably technically correct, it wouldn't help in practice<o:p></o:p></p>
<p>diff --git a/src/projections/aitoff.cpp b/src/projections/aitoff.cpp<br>
index 9d060999f..cce319157 100644<br>
--- a/src/projections/aitoff.cpp<br>
+++ b/src/projections/aitoff.cpp<br>
@@ -59,6 +59,11 @@ static PJ_XY aitoff_s_forward(PJ_LP lp, PJ *P) { /* Spheroidal, forward */<br>
struct pj_opaque *Q = static_cast<struct pj_opaque *>(P->opaque);<br>
double c, d;<br>
<br>
+ if( fabs(lp.lam) > 2 * M_PI ) {<br>
+ proj_errno_set(<br>
+ P, PROJ_ERR_COORD_TRANSFM_OUTSIDE_PROJECTION_DOMAIN);<br>
+ return xy;<br>
+ }<br>
c = 0.5 * lp.lam;<br>
d = acos(cos(lp.phi) * cos(c));<br>
if (d != 0.0) { /* basic Aitoff */<o:p></o:p></p>
<p>The only solution I can offer is that you must explictly clip your source geometry to the validity range of wintri with -clipsrc -360 -90 360 90.<o:p></o:p></p>
<p>Even<o:p></o:p></p>
<p><o:p> </o:p></p>
<div>
<p class="MsoNormal">Le 30/04/2023 à 03:11, Joaquim Manuel Freire Luís via gdal-dev a écrit :<o:p></o:p></p>
</div>
<blockquote style="margin-top:5.0pt;margin-bottom:5.0pt">
<p class="MsoNormal">Even, this is a kind of continuation of the subject that I brought up in
<a href="https://github.com/OSGeo/gdal/issues/7644">https://github.com/OSGeo/gdal/issues/7644</a><o:p></o:p></p>
<p class="MsoNormal"> <o:p></o:p></p>
<p class="MsoNormal">The file size limitation of attachments doesn’t help to make this case easily reproducible, but the story is, I can make that rectangular map using gdalwarp and ogr2ogr that use the +over option and +proj=vandg<o:p></o:p></p>
<p class="MsoNormal">(see <a href="https://forum.generic-mapping-tools.org/t/best-projection-for-rectangular-world-map/3715/82?u=joaquim">
https://forum.generic-mapping-tools.org/t/best-projection-for-rectangular-world-map/3715/82?u=joaquim</a>)<o:p></o:p></p>
<p class="MsoNormal"> <o:p></o:p></p>
<p class="MsoNormal">But when I try the same with +proj=wintry the gdalwarp op works well but the ogr2ogr doesn’t.<o:p></o:p></p>
<p class="MsoNormal">(see <a href="https://forum.generic-mapping-tools.org/t/best-projection-for-rectangular-world-map/3715/85?u=joaquim">
https://forum.generic-mapping-tools.org/t/best-projection-for-rectangular-world-map/3715/85?u=joaquim</a>)<o:p></o:p></p>
<p class="MsoNormal"> <o:p></o:p></p>
<p class="MsoNormal">You can more less reproduce this case using this file<o:p></o:p></p>
<p class="MsoNormal"><a href="http://fct-gmt.ualg.pt/tmp/cl540.gpkg">http://fct-gmt.ualg.pt/tmp/cl540.gpkg</a><o:p></o:p></p>
<p class="MsoNormal"> <o:p></o:p></p>
<p class="MsoNormal">ogr2ogr -t_srs "+proj=wintri +over" cl540_wintri.gpkg cl540.gpkg<o:p></o:p></p>
<p class="MsoNormal"> <o:p></o:p></p>
<p class="MsoNormal">Why would that be in this case?<o:p></o:p></p>
<p class="MsoNormal"> <o:p></o:p></p>
<p class="MsoNormal">Joaquim<o:p></o:p></p>
<p class="MsoNormal"><span style="mso-ligatures:none"><br>
<br>
<o:p></o:p></span></p>
<pre>_______________________________________________<o:p></o:p></pre>
<pre>gdal-dev mailing list<o:p></o:p></pre>
<pre><a href="mailto:gdal-dev@lists.osgeo.org">gdal-dev@lists.osgeo.org</a><o:p></o:p></pre>
<pre><a href="https://lists.osgeo.org/mailman/listinfo/gdal-dev">https://lists.osgeo.org/mailman/listinfo/gdal-dev</a><o:p></o:p></pre>
</blockquote>
<pre>-- <o:p></o:p></pre>
<pre><a href="http://www.spatialys.com">http://www.spatialys.com</a><o:p></o:p></pre>
<pre>My software is free, but my time generally not.<o:p></o:p></pre>
</div>
</body>
</html>