Implement OGC (edge of pixel) extent calculation

szekerest · szekerest · commit 4701966fd20f · 2025-11-05T22:32:16.000+01:00
diff --git a/src/hittest.c b/src/hittest.c
@@ -187,8 +187,8 @@ int msHitTestLayer(mapObj *map, layerObj *layer, layer_hittest *hittest) {
             &searchrect); /* project the searchrect to source coords */
     } else {
       searchrect.minx = searchrect.miny = 0;
-      searchrect.maxx = map->width - 1;
-      searchrect.maxy = map->height - 1;
+      searchrect.maxx = map->width;
+      searchrect.maxy = map->height;
     }
 #ifdef USE_GEOS
     msInitShape(&searchpoly);
diff --git a/src/interpolation.c b/src/interpolation.c
@@ -138,8 +138,8 @@ int msInterpolationDataset(mapObj *map, imageObj *image,
     }
   } else {
     searchrect.minx = searchrect.miny = 0;
-    searchrect.maxx = map->width - 1;
-    searchrect.maxy = map->height - 1;
+    searchrect.maxx = map->width;
+    searchrect.maxy = map->height;
   }
 
   layer->project =
diff --git a/src/mapchart.c b/src/mapchart.c
@@ -688,8 +688,8 @@ int msDrawChartLayer(mapObj *map, layerObj *layer, imageObj *image) {
       searchrect = map->extent;
     else {
       searchrect.minx = searchrect.miny = 0;
-      searchrect.maxx = map->width - 1;
-      searchrect.maxy = map->height - 1;
+      searchrect.maxx = map->width;
+      searchrect.maxy = map->height;
     }
 
     if ((map->projection.numargs > 0) && (layer->projection.numargs > 0))
diff --git a/src/mapcluster.c b/src/mapcluster.c
@@ -1082,8 +1082,8 @@ int RebuildClusters(layerObj *layer, int isQuery) {
     searchrect = map->extent;
   else {
     searchrect.minx = searchrect.miny = 0;
-    searchrect.maxx = map->width - 1;
-    searchrect.maxy = map->height - 1;
+    searchrect.maxx = map->width;
+    searchrect.maxy = map->height;
   }
 
   if (searchrect.minx == layerinfo->searchRect.minx &&
diff --git a/src/mapdraw.c b/src/mapdraw.c
@@ -225,6 +225,11 @@ imageObj *msPrepareImage(mapObj *map, int allow_nonsquare) {
   } else
     map->cellsize = msAdjustExtent(&(map->extent), map->width, map->height);
 
+  if (map->cellsize == 0) {
+    msSetError(MS_MISCERR, "Invalid map extent or image size.", "msAdjustExtent()");
+    return(NULL);
+  }
+
   status = msCalculateScale(map->extent, map->units, map->width, map->height,
                             map->resolution, &map->scaledenom);
   if (status != MS_SUCCESS) {
@@ -1173,8 +1178,8 @@ int msDrawVectorLayer(mapObj *map, layerObj *layer, imageObj *image) {
 
   } else {
     searchrect.minx = searchrect.miny = 0;
-    searchrect.maxx = map->width - 1;
-    searchrect.maxy = map->height - 1;
+    searchrect.maxx = map->width;
+    searchrect.maxy = map->height;
   }
 
   status = msLayerWhichShapes(layer, searchrect, MS_FALSE);
diff --git a/src/mapgraticule.c b/src/mapgraticule.c
@@ -663,8 +663,8 @@ msGraticuleLayerGetIntersectionPoints(mapObj *map, layerObj *layer) {
     searchrect = map->extent;
   else {
     searchrect.minx = searchrect.miny = 0;
-    searchrect.maxx = map->width - 1;
-    searchrect.maxy = map->height - 1;
+    searchrect.maxx = map->width;
+    searchrect.maxy = map->height;
   }
 
   if ((map->projection.numargs > 0) && (layer->projection.numargs > 0))
diff --git a/src/mapobject.c b/src/mapobject.c
@@ -318,7 +318,7 @@ int msMapScaleExtent(mapObj *map, double zoomfactor, double minscaledenom,
 
   if (minscaledenom > 0 || maxscaledenom > 0) {
     /* ensure we are within the valid scale domain */
-    md = (map->width - 1) /
+    md = (map->width) /
          (map->resolution * msInchesPerUnit(map->units, center_y));
     if (minscaledenom > 0 && geo_width < minscaledenom * md)
       geo_width = minscaledenom * md;
@@ -422,14 +422,14 @@ int msMapComputeGeotransform(mapObj *map)
 
 {
   /* Do we have all required parameters? */
-  if (map->extent.minx == map->extent.maxx || map->width <= 1 ||
-      map->height <= 1)
+  if (map->extent.minx == map->extent.maxx || map->width <= 0 ||
+      map->height <= 0)
     return MS_FAILURE;
 
   const double geo_width = map->extent.maxx - map->extent.minx;
   const double geo_height = map->extent.maxy - map->extent.miny;
-  return msMapComputeGeotransformEx(map, geo_width / (map->width - 1),
-                                    geo_height / (map->height - 1));
+  return msMapComputeGeotransformEx(map, geo_width / map->width,
+                                    geo_height / map->height);
 }
 
 /************************************************************************/
diff --git a/src/mapscale.c b/src/mapscale.c
@@ -104,7 +104,7 @@ int msCalculateScale(rectObj extent, int units, int width, int height,
   case (MS_INCHES):
   case (MS_FEET):
     center_y = (extent.miny + extent.maxy) / 2.0;
-    md = (width - 1) /
+    md = (width) /
          (resolution *
           msInchesPerUnit(
               units, center_y)); /* remember, we use a pixel-center to
@@ -593,7 +593,7 @@ double GetDeltaExtentsUsingScale(double scale, int units, double centerLat,
   case (MS_FEET):
     /* remember, we use a pixel-center to pixel-center extent, hence the width-1
      */
-    md = (width - 1) / (resolution * msInchesPerUnit(units, centerLat));
+    md = (width) / (resolution * msInchesPerUnit(units, centerLat));
     dfDelta = md * scale;
     break;
 
diff --git a/src/maptemplate.c b/src/maptemplate.c
@@ -219,13 +219,13 @@ int setExtent(mapservObj *mapserv) {
         msInchesPerUnit(mapserv->map->units,
                         0); /* user supplied a point and a scale denominator */
     mapserv->map->extent.minx =
-        mapserv->mappnt.x - cellsize * (mapserv->map->width - 1) / 2.0;
+        mapserv->mappnt.x - cellsize * (mapserv->map->width) / 2.0;
     mapserv->map->extent.miny =
-        mapserv->mappnt.y - cellsize * (mapserv->map->height - 1) / 2.0;
+        mapserv->mappnt.y - cellsize * (mapserv->map->height) / 2.0;
     mapserv->map->extent.maxx =
-        mapserv->mappnt.x + cellsize * (mapserv->map->width - 1) / 2.0;
+        mapserv->mappnt.x + cellsize * (mapserv->map->width) / 2.0;
     mapserv->map->extent.maxy =
-        mapserv->mappnt.y + cellsize * (mapserv->map->height - 1) / 2.0;
+        mapserv->mappnt.y + cellsize * (mapserv->map->height) / 2.0;
     break;
   default: /* use the default in the mapfile if it exists */
     if ((mapserv->map->extent.minx == mapserv->map->extent.maxx) &&
diff --git a/src/maptile.c b/src/maptile.c
@@ -345,7 +345,7 @@ int msTileSetExtent(mapservObj *msObj) {
 #ifdef USE_TILE_API
 
   mapObj *map = msObj->map;
-  double dx, dy, buffer;
+  double buffer;
   tileParams params;
 
   /* Read the tile-mode map file parameters */
@@ -468,18 +468,6 @@ int msTileSetExtent(mapservObj *msObj) {
     msDebug("msTileSetExtent(): buffered image size (%d x %d)\n", map->width,
             map->height);
 
-  /*
-  ** Adjust the extents inwards by 1/2 pixel so they are from
-  ** center-of-pixel to center-of-pixel, instead of edge-to-edge.
-  ** This is the way mapserver does it.
-  */
-  dx = (map->extent.maxx - map->extent.minx) / map->width;
-  map->extent.minx += dx * 0.5;
-  map->extent.maxx -= dx * 0.5;
-  dy = (map->extent.maxy - map->extent.miny) / map->height;
-  map->extent.miny += dy * 0.5;
-  map->extent.maxy -= dy * 0.5;
-
   /*
   ** Ensure the labelcache buffer is greater than the tile buffer.
   */
diff --git a/src/maputil.c b/src/maputil.c
@@ -1036,7 +1036,7 @@ int msAdjustImage(rectObj rect, int *width, int *height) {
 double msAdjustExtent(rectObj *rect, int width, int height) {
   double cellsize, ox, oy;
 
-  if (width == 1 || height == 1)
+  if (width == 0 || height == 0)
     return 0;
 
   cellsize = MS_MAX(MS_CELLSIZE(rect->minx, rect->maxx, width),
@@ -1045,9 +1045,9 @@ double msAdjustExtent(rectObj *rect, int width, int height) {
   if (cellsize <= 0) /* avoid division by zero errors */
     return (0);
 
-  ox = MS_MAX(((width - 1) - (rect->maxx - rect->minx) / cellsize) / 2,
+  ox = MS_MAX(((width) - (rect->maxx - rect->minx) / cellsize) / 2,
               0); /* these were width-1 and height-1 */
-  oy = MS_MAX(((height - 1) - (rect->maxy - rect->miny) / cellsize) / 2, 0);
+  oy = MS_MAX(((height) - (rect->maxy - rect->miny) / cellsize) / 2, 0);
 
   rect->minx = rect->minx - ox * cellsize;
   rect->miny = rect->miny - oy * cellsize;
diff --git a/src/mapuvraster.cpp b/src/mapuvraster.cpp
@@ -689,9 +689,9 @@ int msUVRASTERLayerWhichShapes(layerObj *layer, rectObj rect, int isQuery) {
   map_tmp->extent.miny =
       rect.miny - (0.5 * map_cellsize) + (0.5 * map_tmp->cellsize);
   map_tmp->extent.maxx =
-      map_tmp->extent.minx + ((width - 1) * map_tmp->cellsize);
+      map_tmp->extent.minx + (width * map_tmp->cellsize);
   map_tmp->extent.maxy =
-      map_tmp->extent.miny + ((height - 1) * map_tmp->cellsize);
+      map_tmp->extent.miny + (height * map_tmp->cellsize);
 
   if (bHasLonWrap && dfLonWrap == 180.0) {
     if (map_tmp->extent.minx >= 180) {
@@ -863,9 +863,9 @@ int msUVRASTERLayerGetShape(layerObj *layer, shapeObj *shape,
 
   const int x = static_cast<int>(raster_off % uvlinfo->width);
   const int y = static_cast<int>(raster_off / uvlinfo->width);
-  point.x = Pix2Georef(x, 0, uvlinfo->width - 1, uvlinfo->extent.minx,
+  point.x = Pix2Georef(x, 0, uvlinfo->width, uvlinfo->extent.minx,
                        uvlinfo->extent.maxx, MS_FALSE);
-  point.y = Pix2Georef(y, 0, uvlinfo->height - 1, uvlinfo->extent.miny,
+  point.y = Pix2Georef(y, 0, uvlinfo->height, uvlinfo->extent.miny,
                        uvlinfo->extent.maxy, MS_TRUE);
   if (layer->debug == MS_DEBUGLEVEL_VVV)
     msDebug("msUVRASTERLayerWhichShapes(): shapeindex: %ld, x: %g, y: %g\n",
diff --git a/src/mapwcs.cpp b/src/mapwcs.cpp
@@ -2674,9 +2674,9 @@ this request. Check wcs/ows_enable_request settings.",
       /* recompute bounding box so we get exactly the origin and
          resolution requested. */
       params->bbox.maxx =
-          params->bbox.minx + (params->width - 1) * params->resx;
+          params->bbox.minx + params->width * params->resx;
       params->bbox.miny =
-          params->bbox.maxy - (params->height - 1) * params->resy;
+          params->bbox.maxy - params->height * params->resy;
     }
   }
 
diff --git a/src/mapwms.cpp b/src/mapwms.cpp
@@ -1925,25 +1925,6 @@ this request. Check wms/ows_enable_request settings.",
     msFreeCharArray(tokens, n);
   }
 
-  /*
-  ** WMS extents are edge to edge while MapServer extents are center of
-  ** pixel to center of pixel.  Here we try to adjust the WMS extents
-  ** in by half a pixel.  We wait till here because we want to ensure we
-  ** are doing this in terms of the correct WIDTH and HEIGHT.
-  */
-  if (adjust_extent && map->width > 1 && map->height > 1 &&
-      !bbox_pixel_is_point) {
-    double dx, dy;
-
-    dx = (map->extent.maxx - map->extent.minx) / map->width;
-    map->extent.minx += dx * 0.5;
-    map->extent.maxx -= dx * 0.5;
-
-    dy = (map->extent.maxy - map->extent.miny) / map->height;
-    map->extent.miny += dy * 0.5;
-    map->extent.maxy -= dy * 0.5;
-  }
-
   if (request && strcasecmp(request, "DescribeLayer") != 0) {
     if (!srsfound) {
       if (nVersion >= OWS_1_3_0)
diff --git a/src/mapwmslayer.c b/src/mapwmslayer.c
@@ -650,17 +650,8 @@ static int msBuildWMSLayerURL(mapObj *map, layerObj *lp, int nRequestType,
     return MS_FAILURE;
   }
 
-  /* ------------------------------------------------------------------
-   * Adjust for MapServer EXTENT being center of pixel and WMS BBOX being
-   * edge of pixel (#2843).
-   * ------------------------------------------------------------------ */
   bbox = map->extent;
 
-  bbox.minx -= map->cellsize * 0.5;
-  bbox.maxx += map->cellsize * 0.5;
-  bbox.miny -= map->cellsize * 0.5;
-  bbox.maxy += map->cellsize * 0.5;
-
   /* -------------------------------------------------------------------- */
   /*      Reproject if needed.                                            */
   /* -------------------------------------------------------------------- */

Original file line number	Diff line number	Diff line change
`@@ -138,8 +138,8 @@ int msInterpolationDataset(mapObj map, imageObj image,`
`138`	`138`	`}`
`139`	`139`	`} else {`
`140`	`140`	`searchrect.minx = searchrect.miny = 0;`
`141`		`- searchrect.maxx = map->width - 1;`
`142`		`- searchrect.maxy = map->height - 1;`
	`141`	`+ searchrect.maxx = map->width;`
	`142`	`+ searchrect.maxy = map->height;`
`143`	`143`	`}`
`144`	`144`
`145`	`145`	`layer->project =`