UVRaster: support layers with lon_wrap=180 projection (#5502)

mapproject.c

@@ -942,6 +942,32 @@ msProjectRectAsPolygon(projectionObj *in, projectionObj *out,
 #endif
 }
 
+/************************************************************************/
+/*                        msProjectHasLonWrap()                         */
+/************************************************************************/
+
+int msProjectHasLonWrap(projectionObj *in, double* pdfLonWrap)
+{
+    int i;
+    if( pdfLonWrap )
+        *pdfLonWrap = 0;
+#if USE_PROJ
+    if( !pj_is_latlong(in->proj) )
+        return MS_FALSE;
+#endif
+    for( i = 0; i < in->numargs; i++ )
+    {
+        if( strncmp(in->args[i], "lon_wrap=",
+                    strlen("lon_wrap=")) == 0 )
+        {
+            if( pdfLonWrap )
+                *pdfLonWrap = atof(in->args[i] + strlen("lon_wrap="));
+            return MS_TRUE;
+        }
+    }
+    return MS_FALSE;
+}
+
 /************************************************************************/
 /*                           msProjectRect()                            */
 /************************************************************************/
@@ -954,6 +980,7 @@ int msProjectRect(projectionObj *in, projectionObj *out, rectObj *rect)
   char *over = "+over";
   int ret;
   projectionObj in_over,out_over,*inp,*outp;
+  double dfLonWrap = 0.0;
 
 #if USE_PROJ
   /* Detect projecting from north polar stereographic to longlat */
@@ -988,6 +1015,14 @@ int msProjectRect(projectionObj *in, projectionObj *out, rectObj *rect)
   }
 #endif
 
+  if(in && msProjectHasLonWrap(in, &dfLonWrap) && dfLonWrap == 180.0) {
+    inp = in;
+    outp = out;
+    if( rect->maxx > 180.0 ) {
+      rect->minx = -180.0;
+      rect->maxx = 180.0;
+    }
+  }
   /* 
    * Issue #4892: When projecting a rectangle we do not want proj to wrap resulting
    * coordinates around the dateline, as in practice a requested bounding box of
@@ -997,19 +1032,21 @@ int msProjectRect(projectionObj *in, projectionObj *out, rectObj *rect)
    *  To enforce this, we clone the input projections and add the "+over" proj 
    *  parameter in order to disable dateline wrapping.
    */ 
-  if(out) {
-    msInitProjection(&out_over);
-    msCopyProjectionExtended(&out_over,out,&over,1);
-    outp = &out_over;
-  } else {
-    outp = out;
-  }
-  if(in) {
-    msInitProjection(&in_over);
-    msCopyProjectionExtended(&in_over,in,&over,1);
-    inp = &in_over;
-  } else {
-    inp = in;
+  else {
+    if(out) {
+      msInitProjection(&out_over);
+      msCopyProjectionExtended(&out_over,out,&over,1);
+      outp = &out_over;
+    } else {
+      outp = out;
+    }
+    if(in) {
+      msInitProjection(&in_over);
+      msCopyProjectionExtended(&in_over,in,&over,1);
+      inp = &in_over;
+    } else {
+      inp = in;
+    }
   }
   ret = msProjectRectAsPolygon(inp,outp, rect );
   if(in)

mapproject.h

@@ -108,6 +108,8 @@ extern "C" {
 
   /*utility functions */
   MS_DLL_EXPORT int GetMapserverUnitUsingProj(projectionObj *psProj);
+
+  int msProjectHasLonWrap(projectionObj *in, double* pdfLonWrap);
 #endif
 
 #ifdef __cplusplus

mapresample.c

@@ -1136,18 +1136,8 @@ static int msTransformMapToSource( int nDstXSize, int nDstYSize,
   /* -------------------------------------------------------------------- */
   if( bOutInit && pj_is_latlong(psSrcProj->proj) )
   {
-      int bHasLonWrap = MS_FALSE;
       double dfLonWrap = 0;
-      for( i = 0; i < psSrcProj->numargs; i++ )
-      {
-          if( strncmp(psSrcProj->args[i], "lon_wrap=",
-                      strlen("lon_wrap=")) == 0 )
-          {
-              bHasLonWrap = MS_TRUE;
-              dfLonWrap = atof( psSrcProj->args[i] + strlen("lon_wrap=") );
-              break;
-          }
-      }
+      int bHasLonWrap = msProjectHasLonWrap(psSrcProj, &dfLonWrap);
 
       if( bHasLonWrap )
       {

mapuvraster.c

@@ -345,7 +345,8 @@ int msUVRASTERLayerWhichShapes(layerObj *layer, rectObj rect, int isQuery)
   int width, height, u_src_off, v_src_off, i, x, y;
   char   **alteredProcessing = NULL, *saved_layer_mask;
   char **savedProcessing = NULL;
-
+  int bHasLonWrap = MS_FALSE;
+  double dfLonWrap = 0.0;
 
   if (layer->debug)
     msDebug("Entering msUVRASTERLayerWhichShapes().\n");
@@ -389,10 +390,6 @@ int msUVRASTERLayerWhichShapes(layerObj *layer, rectObj rect, int isQuery)
   map_cellsize = MS_MAX(MS_CELLSIZE(rect.minx, rect.maxx,layer->map->width),
                         MS_CELLSIZE(rect.miny,rect.maxy,layer->map->height));
   map_tmp->cellsize = map_cellsize*spacing;
-
-  if (layer->debug)
-    msDebug("msUVRASTERLayerWhichShapes(): width: %d, height: %d, cellsize: %g\n",
-            width, height, map_tmp->cellsize);
 
   /* Initialize our dummy map */
   MS_INIT_COLOR(map_tmp->imagecolor, 255,255,255,255);
@@ -419,7 +416,55 @@ int msUVRASTERLayerWhichShapes(layerObj *layer, rectObj rect, int isQuery)
   map_tmp->extent.maxy = map_tmp->extent.miny+((height-1)*map_tmp->cellsize);
   map_tmp->gt.rotation_angle = 0.0;
 
-   msCopyProjection(&map_tmp->projection, &layer->projection);
+  /* Custom msCopyProjection() that removes lon_wrap parameter */
+  {
+#ifdef USE_PROJ
+    int i;
+
+    map_tmp->projection.numargs = 0;
+    map_tmp->projection.gt = layer->projection.gt;
+    map_tmp->projection.automatic = layer->projection.automatic;
+
+    for (i = 0; i < layer->projection.numargs; i++) {
+      if( strncmp(layer->projection.args[i], "lon_wrap=",
+                  strlen("lon_wrap=")) == 0 ) {
+        bHasLonWrap = MS_TRUE;
+        dfLonWrap = atof( layer->projection.args[i] + strlen("lon_wrap=") );
+      }
+      else {
+        map_tmp->projection.args[map_tmp->projection.numargs ++] =
+            msStrdup(layer->projection.args[i]);
+      }
+    }
+    if (map_tmp->projection.numargs != 0) {
+      msProcessProjection(&(map_tmp->projection));
+    }
+#endif
+    map_tmp->projection.wellknownprojection = layer->projection.wellknownprojection;
+  }
+
+  if( bHasLonWrap && dfLonWrap == 180.0) {
+    if( map_tmp->extent.minx >= 180 ) {
+      /* Request on the right half of the shifted raster (= western hemisphere) */
+      map_tmp->extent.minx -= 360;
+      map_tmp->extent.maxx -= 360;
+    }
+    else if( map_tmp->extent.maxx >= 180.0 ) {
+      /* Request spanning on the 2 hemispheres => drawing whole planet */
+      /* Take only into account vertical resolution, as horizontal one */
+      /* will be unreliable (assuming square pixels...) */
+      map_cellsize = MS_CELLSIZE(rect.miny,rect.maxy,layer->map->height);
+      map_tmp->cellsize = map_cellsize*spacing;
+
+      width = 360.0 / map_tmp->cellsize;
+      map_tmp->extent.minx = -180.0+(0.5*map_tmp->cellsize);
+      map_tmp->extent.maxx = 180.0-(0.5*map_tmp->cellsize);
+    }
+  }
+
+  if (layer->debug)
+    msDebug("msUVRASTERLayerWhichShapes(): width: %d, height: %d, cellsize: %g\n",
+            width, height, map_tmp->cellsize);
 
   if (layer->debug == 5)
     msDebug("msUVRASTERLayerWhichShapes(): extent: %g %g %g %g\n",

msautotest/wxs/wms_uvraster_lonwrap180.map

@@ -0,0 +1,74 @@
+#
+# Test WMS
+#
+# REQUIRES: INPUT=GDAL OUTPUT=PNG SUPPORTS=WMS
+# 
+# RUN_PARMS: wms_uvraster_lonwrap180_whole_world.png [MAPSERV] QUERY_STRING="map=[MAPFILE]&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=-90,-180,90,180&CRS=EPSG:4326&WIDTH=744&HEIGHT=372&LAYERS=test&STYLES=&FORMAT=image/png&DPI=120&MAP_RESOLUTION=120&TRANSPARENT=TRUE" > [RESULT_DEMIME]
+#
+# RUN_PARMS: wms_uvraster_lonwrap180_west_hemisphere.png [MAPSERV] QUERY_STRING="map=[MAPFILE]&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=-90,-180,90,0&CRS=EPSG:4326&WIDTH=372&HEIGHT=372&LAYERS=test&STYLES=&FORMAT=image/png&DPI=120&MAP_RESOLUTION=120&TRANSPARENT=TRUE" > [RESULT_DEMIME]
+#
+# RUN_PARMS: wms_uvraster_lonwrap180_east_hemisphere.png [MAPSERV] QUERY_STRING="map=[MAPFILE]&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=-90,0,90,180&CRS=EPSG:4326&WIDTH=372&HEIGHT=372&LAYERS=test&STYLES=&FORMAT=image/png&DPI=120&MAP_RESOLUTION=120&TRANSPARENT=TRUE" > [RESULT_DEMIME]
+#
+# RUN_PARMS: wms_uvraster_lonwrap180_accross_hemisphere.png [MAPSERV] QUERY_STRING="map=[MAPFILE]&SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=-90,-90,90,90&CRS=EPSG:4326&WIDTH=372&HEIGHT=372&LAYERS=test&STYLES=&FORMAT=image/png&DPI=120&MAP_RESOLUTION=120&TRANSPARENT=TRUE" > [RESULT_DEMIME]
+
+MAP
+ NAME uvraster_lonwrap
+ EXTENT -180 -90 180 90
+ MAXSIZE 4096
+ SIZE 500 300
+ SHAPEPATH ./data
+
+ PROJECTION
+  "init=epsg:4326"
+ END
+
+ SYMBOL
+  NAME "arrow"
+  FILLED TRUE
+  POINTS
+   1.2 0
+   1.2 6
+   3 6
+   0 10
+   -3 6
+   -1.2 6
+   -1.2 0
+   1.2 0
+  END
+ END
+
+ WEB
+  METADATA
+    "ows_enable_request" "*"
+  END
+ END
+
+ LAYER
+  DATA "uvraster_lonwrap180.tif"
+  NAME "test"
+  TYPE POINT
+  CONNECTIONTYPE uvraster
+
+  EXTENT 0 -80.1 360 89.9
+
+  PROCESSING "BANDS=1,2"
+  PROCESSING "UV_SPACING=17"
+
+  PROJECTION
+   "proj=longlat"
+   "a=6371229"
+   "b=6371229"
+   "lon_wrap=180"
+  END
+
+  CLASS
+   STYLE
+    SYMBOL "arrow"
+    COLOR 255 0 0
+    ANGLE [uv_angle]
+    SIZE 10
+   END
+  END
+
+ END
+END