Computed CRS EXTENT for EPSG 32661 UPS North - ignore axis orientation - X/Y inversion issue
Hi,
I am trying to upgrade from an older WMS service running mapserver 7.6.4 to a newer machine running mapserver 8.5 (full version below [1] [2])
-
Old service:
mapserv -v[1] MapServer version 7.6.4 OUTPUT=PNG OUTPUT=JPEG OUTPUT=KML SUPPORTS=PROJ SUPPORTS=AGG SUPPORTS=FREETYPE SUPPORTS=CAIRO SUPPORTS=SVG_SYMBOLS SUPPORTS=RSVG SUPPORTS=ICONV SUPPORTS=FRIBIDI SUPPORTS=WMS_SERVER SUPPORTS=WMS_CLIENT SUPPORTS=WFS_SERVER SUPPORTS=WFS_CLIENT SUPPORTS=WCS_SERVER SUPPORTS=SOS_SERVER SUPPORTS=FASTCGI SUPPORTS=THREADS SUPPORTS=GEOS SUPPORTS=POINT_Z_M SUPPORTS=PBF INPUT=JPEG INPUT=POSTGIS INPUT=OGR INPUT=GDAL INPUT=SHAPEFILE -
New service:
mapserv -v[2] MapServer version 8.5-dev PROJ version 9.4 GDAL version 3.10 OUTPUT=PNG OUTPUT=JPEG SUPPORTS=PROJ SUPPORTS=AGG SUPPORTS=FREETYPE SUPPORTS=CAIRO SUPPORTS=SVG_SYMBOLS SUPPORTS=RSVG SUPPORTS=ICONV SUPPORTS=XMP SUPPORTS=FRIBIDI SUPPORTS=WMS_SERVER SUPPORTS=WMS_CLIENT SUPPORTS=WFS_SERVER SUPPORTS=WFS_CLIENT SUPPORTS=WCS_SERVER SUPPORTS=OGCAPI_SERVER SUPPORTS=THREADS SUPPORTS=GEOS SUPPORTS=PBF INPUT=JPEG INPUT=POSTGIS INPUT=OGR INPUT=GDAL INPUT=SHAPEFILE INPUT=FLATGEOBUF
The dataset used to test the two services are exactly the same and so are the mapfiles (except for path to fonts, symbols, etc)
I reduced the mapfile used on the NEW service, which is equivalent to the OLD service, to a single layer (GTiff in EPSG:32632), and looks like:
MAP
EXTENT 5.38027653341219 60.0865467027648 7.31486520966422 61.536882894348
# extent 5.26749797619507 60.2993491041863 7.36602770443023 61.3239944400882
IMAGETYPE "png"
NAME "NBS"
SIZE 800 600
STATUS ON
UNITS METERS
OUTPUTFORMAT
NAME "png"
MIMETYPE "image/png"
DRIVER "AGG/PNG"
EXTENSION "png"
IMAGEMODE RGB
TRANSPARENT FALSE
END # OUTPUTFORMAT
PROJECTION
"init=EPSG:4326"
END # PROJECTION
LEGEND
KEYSIZE 20 10
KEYSPACING 5 5
LABEL
SIZE 10
OFFSET 0 0
SHADOWSIZE 1 1
END # LABEL
STATUS OFF
END # LEGEND
QUERYMAP
STATUS OFF
STYLE HILITE
END # QUERYMAP
SCALEBAR
INTERVALS 4
LABEL
SIZE 10
OFFSET 0 0
SHADOWSIZE 1 1
END # LABEL
SIZE 200 3
STATUS OFF
UNITS MILES
END # SCALEBAR
WEB
METADATA
"wms_srs" "EPSG:32632 EPSG:4326 EPSG:32661"
"wms_onlineresource" "https://api.get_wms_prod/f521c540-486f-42f9-9871-df3a54bf172c/wms"
"wms_enable_request" "*"
"wms_abstract" "Each of the satellites in the SENTINEL-2 mission carries a single payload: the Multi-Spectral Instrument (MSI)."
"wms_bbox_extended" "true"
"wms_title" "S2A_MSIL1C_20250120T111401_N0511_R137_T32VLN_20250120T112920"
END # METADATA
END # WEB
LAYER
DATA "false_color_glacier.tiff"
PROJECTION
"init=epsg:32632"
END
EXTENT 299995 6690245 409795 6800045
METADATA
"wms_srs" "EPSG:32632 EPSG:4326 EPSG:32661"
"wms_abstract" "RGB Composite of bands B11, B8, B4"
"wms_layer_group" "/composite"
"wms_extent" "299995 6690245 409795 6800045"
"wms_group_title" "Composite"
"wms_title" "False Color Glacier Composite"
END # METADATA
NAME "false_color_glacier"
STATUS ON
TILEITEM "location"
TYPE RASTER
UNITS METERS
END # LAYER
END # MAP
I can get the NEW service working with all the projection except for the EPSG:32661 (UPS North, which has inverted axis) when I load such layer in any client it returns blank tiles for the area where it is supposed to be - and even in QGIS, loading the layer results in a map that gets distorsion and disappears as soon as I zonn/pan on the canvas.
Trying to debug the issue, and exploring the getcapabilities of the two service, I noticed that the CRS BBOX Extent is OFF on the new service, X/Y are swapped for the projection I need to use: EPSG 32661 (which indeed as inverted axis but seems that my mapserver fails to detect that)
Of course the OLD service has no problem in let any client to render request using EPSG:32661 - see the details inb the relative getCapabilities:
<BoundingBox CRS="EPSG:32661" minx="2.2983e+06" miny="-1.35614e+06" maxx="2.43211e+06" maxy="-1.22239e+06"/>
<BoundingBox CRS="EPSG:32661" minx="-1.35552e+06" miny="2.29873e+06" maxx="-1.22275e+06" maxy="2.43151e+06"/>
Can You plerase help understanding whatr I am missing on the new service to get this to work? AS tesrting environemt on the NEW service I am using mapserver docker image from camptocamp - I built it from mapserver GIT master branch.
I tried some more debug, I puilled an older mapserver docker image: camptocamp/mapserver:7.6-gdal3.8 and rewrote the mapfile as following:
MAP
EXTENT 5.38027653341219 60.0865467027648 7.31486520966422 61.536882894348
# extent 5.26749797619507 60.2993491041863 7.36602770443023 61.3239944400882
IMAGETYPE "png"
NAME "NBS"
# SIZE 800 600
SIZE 650 600
MAXSIZE 4096
STATUS ON
UNITS METERS
CONFIG "PROJ_LIB" "/usr/local/share/proj/"
CONFIG "PROJ_DATA" "/usr/local/share/proj/"
OUTPUTFORMAT
NAME "png"
MIMETYPE "image/png"
DRIVER "AGG/PNG"
EXTENSION "png"
IMAGEMODE RGB
TRANSPARENT FALSE
END # OUTPUTFORMAT
PROJECTION
"init=epsg:4326"
END # PROJECTION
LEGEND
KEYSIZE 20 10
KEYSPACING 5 5
LABEL
SIZE 10
OFFSET 0 0
SHADOWSIZE 1 1
END # LABEL
STATUS OFF
END # LEGEND
QUERYMAP
STATUS OFF
STYLE HILITE
END # QUERYMAP
SCALEBAR
INTERVALS 4
LABEL
SIZE 10
OFFSET 0 0
SHADOWSIZE 1 1
END # LABEL
SIZE 200 3
STATUS OFF
UNITS MILES
END # SCALEBAR
WEB
METADATA
"wms_srs" "EPSG:32632 EPSG:4326 EPSG:32661"
"wms_onlineresource" "https://mapserver7.wps.met.no/?map=/etc/mapserver/wms_ref.map"
"wms_enable_request" "*"
"wms_abstract" "Each of the satellites in the SENTINEL-2 mission carries a single payload: the Multi-Spectral Instrument (MSI)."
"wms_bbox_extended" "true"
"wms_title" "S2A_MSIL1C_20250120T111401_N0511_R137_T32VLN_20250120T112920"
END # METADATA
END # WEB
LAYER
DATA "false_color_glacier.tiff"
PROJECTION
"init=epsg:32632"
END
EXTENT 299995 6690245 409795 6800045
METADATA
"wms_srs" "EPSG:32632 EPSG:4326 EPSG:32661"
"wms_abstract" "RGB Composite of bands B11, B8, B4"
"wms_layer_group" "/composite"
"wms_extent" "299995 6690245 409795 6800045"
"wms_group_title" "Composite"
"wms_title" "False Color Glacier Composite"
END # METADATA
NAME "false_color_glacier"
STATUS ON
TILEITEM "location"
TYPE RASTER
UNITS METERS
END # LAYER
END # MAP
I added:
CONFIG "PROJ_LIB" "/usr/share/proj/"
CONFIG "PROJ_DATA" "/usr/share/proj/"
and renamed the projection string to lower case:
- at map level:
PROJECTION
"init=epsg:4326"
END # PROJECTION
- at layer level:
PROJECTION
"init=epsg:32632"
END
the getcapabilities for this version+mapfile returns the expected result.
<BoundingBox CRS="EPSG:32661" minx="2.30231e+06" miny="-1.38065e+06" maxx="2.43234e+06" maxy="-1.19785e+06"/>
can you help me undertanding what I am doing wrong with mapserver 8.x ? Am I using a misconfigured mapfile or missing bits in the mapserver configuration?
Thanks!!!
EDIT: I updated to a more recent version of ubuntu/gdal while keeping mapserver to 7.6 (same issue)
EPSG:32661 (Arctic) is marked as swapaxis but not neu in proj4 definition since a while back (unsure exactly what version of proj it did change)
also a fix for qgis have been done for the same ting
https://github.com/qgis/QGIS/pull/60597
I have updated the service running mapserver 7.6 to include python mapscript (which I use to generate and dispatch the mapobject in my use case)
https://mapserver7py.wps.met.no/?map=/etc/mapserver/wms_ref.map&service=WMS&request=GetCapabilities
Any thought on how to proceed? As at the moment at work they require to use 32661 and I am now forced to port my app back to mapserver7.x if the issue is not addressed upstream in Mapserver
@kungzlatan should I try a more recent build of QGIS? I think the issue is on mapserver side as the older version 7.x is rendering correctly - I have been trying also client like OpenLayers - It has the same problem.
Confirmed to be an issue/regression introduced with mapserver 8.x
From the PROJ developers,
I can't find any information abut axis orientation
really ;-) ?
CS[Cartesian,2], AXIS["northing (N)",south, MERIDIAN[180, ANGLEUNIT["degree",0.0174532925199433]], ORDER[1], LENGTHUNIT["metre",1]], AXIS["easting (E)",south, MERIDIAN[90, ANGLEUNIT["degree",0.0174532925199433]], ORDER[2], LENGTHUNIT["metre",1]],This is a northing, easting ordered CRS. This is a MapServer issue not a PROJ one.
@epifanio
This change that was introduced in 8.2. https://github.com/MapServer/MapServer/pull/7096. So the issue was that in earlier versions of MapServer, EPSG:32661 was wrongly treated as having axis=enu.
Regarding the GetCapabilities, I think it is correct as it is in the example for version 8.5 that you provided. According to the WMS specification, the x in <BoundingBox> corresponds to the first axis in the CRS (i.e. y in 32661). See discussion at: https://github.com/MapServer/MapServer/issues/6978
So if we transform the minx and miny coordinates from the BBOX in you example, I get values close to the extent defined in the MAP file:
<BoundingBox CRS="EPSG:32661" minx="-1.35552e+06" miny="2.29873e+06"
minx=-1.35552e+06
miny=2.29873e+06
x=$miny # x in <BoundingBox> corresponds to y in 32661
y=$minx # y in <BoundingBox> corresponds to x in 32661
echo $x $y | gdaltransform -s_srs EPSG:32661 -t_srs EPSG:4326 # 5.08742708186411 60.3127167369556 0
The GetCapabilities links do not work at the moment, but next time it would be good to write the request with all parameters, expecially including &VERSION=1.3.0 because the axis order thing was changed between versions 1.1.1 and 1.3.0 (now somebody will say that it was not changed but it was only understood in a wrong way before - it means the same from users perspective).
Maybe Mapserver is using &version=1.3.0 by default but it does not hurt to be explicit.
yes, I was forcing the version to be 1.3 - I create the getcapabilities from inside FastAPI via mapscript.
As suggested in https://github.com/MapServer/MapServer/issues/7273 - I rebuilt QGIS using git master - and it is now loading the WMS layer directly from the endpoint even when specifying the CRS EPSG:32661
This is great, as it seems that now both QGIS / Mapserver / GDAL / Proj are aligned in using the same axis convention when dealing with inverted axis projection.
As mentioned in https://github.com/MapServer/MapServer/issues/7273 the problem will now arise when loading the WMS into a JS client - to see if it can handle EPSG:32661 correctly.
I have some OL code to work wit proj-js and switch between EPSG:4326 and EPSG:32661 but this was using mapserver7.x - the same client code fails to render the WMS when provided by mapserver V 8.x
I will try to make sure I have the latest version of both OL and proj-js and try again - and report back here @PatrikSylve @jratike80 - thanks a lot for looking itno this :)
this is now getting interesting, I changed the FastAPI mapserver API to use a jinja template and create a mapfile (based on mapserver 7.x ) which is then served via mapscript - this is using mapscript + mapserver V7.x
A temporary endpoint: https://map.wps.met.no/get_wms/26c3a197-c590-4e2b-bd44-ff213b56a1ce/wms
this works as expected in the "old" QGIS - I tried V3.34.4-Prizren - but not in QGIS master (I guess this was gonna be expected)
I guess my only alternative to have a stable service that will work with any client is to implement some logic to detect the client that is sending the request ..
I may try to build 2 jinja templates (mapserver V7 + V8) and instead of dispatching the service via FastAPI+mapscript I will have to redirect the response to a mapserver fastcgi service. like having 2 mapserver fcgi services running in different containers (one for V7 and an other for V8) and let the FastAPI redirect to the mapserv version based on some information from the client request :( it seems painful :(
Hi,
I have QGIS 3.43.0, installed with OSGeo4W. I set the QGIS project CRS into EPSG:32661 and asked WMS to use the same CRS. Zoom to layer works OK but zooming in returns results that look unexpected to me and they do not fill the bbox that I zoomed in. But maybe it is normal because the uncommon asis orientation of this coordinate system https://geohubkenya.wordpress.com/2019/08/15/ups-coordinate-system/
One request captured from QGIS as an example.
https://map.wps.met.no/get_wms/26c3a197-c590-4e2b-bd44-ff213b56a1ce/wms?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&BBOX=1120010%2C-692280%2C1702890%2C90152.19999999999709&CRS=EPSG%3A32661&WIDTH=159&HEIGHT=120&LAYERS=composite&STYLES=&FORMAT=image%2Fpng&DPI=144&MAP_RESOLUTION=144&FORMAT_OPTIONS=dpi%3A144&TRANSPARENT=TRUE