Axis order in WMS 1.3.0

[elzouavo]

WMS 1.3.0 standard says "Coordinates shall be listed in the order defined by the CRS and shall be mapped appropriately...". MapServer 5.6.5 hard codes that EPSG codes between 4000 and 5000 all have axes with an inverted order.

Currently, I have to update my data servers to add the EPSG codes for French overseas territories (with GRS80). In addition I have to add the codes for the INSPIRE directive (for France metropolitan).

It turns out that EPSG codes:4559 (French Antilles) and EPSG:4471 (Mayotte) are in the range 4000-5000 but have not reversed axes order in the EPSG database (which is correct). So, customers who question my WMS servers must make wrong requests to have a right result!

With EPSG codes recommended by the European INSPIRE directive, only the EPSG:4258 is properly handled by MapServer 5.6.5 (range 4000-5000). All other concerning France (EPSG:3042, EPSG:3043 and EPSG:3044) and Pan-Europe EPSG:3034 and EPSG:3035, by definition, have reverse axes order. However, they are not in the range 4000-5000, MapServer takes as coordinate systems with "normal" axes order (longitude, latitude). Again, customers who request my WMS servers must do (if possible ...) false queries for correct answers!

I realize that WMS 1.3.0 standard (and WFS 1.1.0) brings news in coordinate systems which are not easy to integrate. OpenLayers 2.10 has added parameter yx for "layers" (which allows great flexibility).

I'm not a developer so I have no advice to give. But would not it be interesting to add a parameter in the configuration of MapServer to point to a list of codes that must be considered with inverse axes order. Example: CONFIG "CRS_YX" "C:/somedir/crs_yx.txt"

Best regards.

[rouault]

An idea to implement this correctly : when mapserver is compiled with GDAL/OGR support, it could use the OSREPSGTreatsAsLatLong() function to know more reliably the axis order taken from the EPSG database. The only drawback is potential longer time execution as you need to instantiate the OGRSpatialReference object from the the gdal/data/*.csv files derived from the EPSG database.

[assefa]

I was not aware of these issues.

I am not really familiar with the projection details so some questions: would OSREPSGTreatsAsLatLong return true for codes that use inverted axis? Is the axis ordering info read from the espg file?

[rouault]

Yes and yes. Basically you do

OSRImportFromEPSGA(hSRS, nEPSGCode); see ​http://gdal.org/ogr/classOGRSpatialReference.html#aa6965a1df98cdc673dfb20697eab613 int bIsInvertedAxisOrder = OSREPSGTreatsAsLatLong(hSRS); ​http://gdal.org/ogr/classOGRSpatialReference.html#519cce70fd39f974c61342bf73ab01ad

For OSREPSGTreatsAsLatLong() apparently it only deals with geographic systems and not projected systems. So perhaps the logic should be customized to deal with that case.

[elzouavo]

I think the implementation would provide maximum flexibility and not include everything in the compilation. MapServer can be used as a server WMS (if we limit to WMS) as well as a client. And, the MapServer WMS client should be able to adapt to query a WMS server 1.3.0 which does not properly deals with reversed axis order. As OpenLayers 2.10 must do today to ask MapServer 5.6.5 for many CRS. Actually it is not only geographic coordinate systems that have reversed axes order. Most INSPIRE systems I mentioned are projected systems.

[assefa]

Sorry for the delay. I am working on this as part of the 6.0 release. I will port any modifications to the 5.6.x branch.

Checking this earlier, I was able to test for wms cases the use of OSRxxx functions mentioned in this bug. It works as excepted, provided that all the utility files are found (gss.csv, ....) As mentionned alos in the bug, the use of the function OSREPSGTreatsAsLatLong will only work for non projected coordinates. Is there an equivalent function that deal with projected system? I am not really sure I understand why that function need to make sure that It is a non-projected coordinates. I was candidly assuming that I would be able to define my axis order somehow in coordinate_axis.csv for all the projections that I wanted "inverted"

[rouault]

No, AFAIK, there's no equivalent of OSREPSGTreatsAsLatLong() for projected CRS in OSR. Indeed, the name of the function makes it difficult to extend its behaviour. Perhaps a OSREPSGIsInvertedAxisOrder() (not sure the naming is really good) would be needed there. Perhaps Frank has thoughts on this.

[warmerdam]

OGRSpatialReference() also has the GetAxis() method which can be used to identify the more general case of axis handling for projected coordinate systems. The EPSGTreatsAsLatLong() method was mostly intended as a stopgap.

I'm pretty nervous about MapServer having to make extensive use of OGRSpatialReference for common cases of EPSG handling. Perhaps it would be better to pre-generate a small axis-orientation table and embed this in MapServer along with a script to regenerate it at will?

[assefa]

ok I will try it with a static table in a .h file that can be included. I think the table will only need the epsg code and a 1 or 0 for inverted/not inverted. If the epsg code is not found, we will pass through the current mechanism that epsg between 4000 and 5000 are inverted. I will initially try it with wms.

[assefa]

added a static table to define some espg codes and their axis order. r11155. There is a mapaxisorder.csv and a mapaxisorder.sh that can be used to generate the mapaxisorder.h that contains the static table. if the espg code is not in this table, the same logic as now apples (espg code between 4000 and 500 are considered to have inverted axis. hopefully this is enough for the cases described here.

Frank do you see any problem using the msIsAxisInverted function in msLoadProjectionString (mapfile.c)?

[warmerdam]

Assefa,

This looks great, thanks!

[elzouavo]

This is perhaps great when MapServer is used as WMS Server, but it seems to be not so great when MapServer is used as WMS Client (to fit the client). How can MapServer fit the WMS client with this implementation ? Thanks for all

[assefa]

closing this. There is a specfic bug for wMS client upgrade to 1.30 (#3039)

[mko]

1 mapaxisorder.(h|csv) describes epsg:4559 having inverted axis order: {4559,1}, , but according to ​http://www.epsg-registry.org/ epsg:4559 has east/north orientation, so the inverted flag should be removed for this code.

---

2 all occurrences of if( atoi(code) >= 4000 && atoi(code) < 5000 ) in mapfile.c should be changed to if( msIsAxisInverted(atoi(code)) )

---

3 MapServer 6 is still not in line with WMS 1.3.0 regarding the axis order. The list of inverted epsg codes is unsatisfying and will result in servers running MapServer 6 correctly handling the axis order (user-compiled version with a fully adjusted mapaxisorder.h), partly correct axis order (some adjustments), the default mapaxisorder.h, different pre-compiled package versions or even worse. Quickly querying the epsg database

SELECT * FROM [Coordinate Reference System] 
WHERE COORD_SYS_CODE IN 
( SELECT CAA.COORD_SYS_CODE FROM [Coordinate Axis] AS CAA 
  INNER JOIN [Coordinate Axis] AS CAB 
  ON CAA.COORD_SYS_CODE = CAB.COORD_SYS_CODE 
  WHERE caa.ORDER=1 AND cab.ORDER=2
  AND ( 
   ( left(CAA.COORD_AXIS_ORIENTATION,5)="north" 
     and left(CAB.COORD_AXIS_ORIENTATION,4)="east" ) 
   or ( left(CAA.COORD_AXIS_ORIENTATION,5)="south" 
     and left(CAB.COORD_AXIS_ORIENTATION,4)="west" ) 
  ) 
);

there should be about 1700 epsg codes with axis order north/east; unaudited mapaxisorder.h attached to this ticket.

In order to make MapServer serving WMS 1.3.0 correctly, we should examine the list of epsg codes and implement them. Otherwise users will have ongoing troubles with WMS 1.3.0.

[mko]

unchecked list of about 1700 epsg codes probably having north/east axes orientation

[assefa]

2 infos coming from the list

• Martin, can you attach the mapaxisorder.csv if you have it
• I am not sure how to confirm this. But I am willing to take the list as is.

Should this be applied this for 6.0, rc2 is tomorrow? I am questioning this mainly because the mechanism used to check this list goes through the table 1 by 1 and is not optimized for a large table. I don't thing it is really costly compared to drawing time though. Also the list is not confirmed although It is the best list available (at least for me since I lack specific knowledge on this matter). Maybe in a point release would be more appropriate.

[assefa]

patch applied in r11646. Thanks

[assefa]

saw this on IRC logs from dmorissette.:

"anther possible optimization (for a future release) could be to use a bit array instead of a table of epsg codes as we have now. EPSG codes must be in the range 0-32767 AFAIK (codes above 32767 are not official EPSG AFAIK)... so that would mean a bit array of 4096 bytes (smaller than the current array of codes) and a lookup would be very fast"

[mko]

What about removing the code >= 4000 && code < 5000 part and searching through the list of codes with swaped axes, only?

[dmorissette]

I don't understand the last comment. I believe Assefa's r11646 (in trunk, to be released in 6.0) already got rid of (code >= 4000 && code < 5000) and uses the list of swapped axes only. Please provide more specific info if there is still an issue with the latest trunk version.

[assefa]

Last comment is not an issue. Martin re-upadted the same patch file that was initially there and indicated that this new patch file was suppressing warnings. I did not have the warnings mentioned locally so I did not apply it. I have communicated on this with Martin.

The comment was valid the first time the patch was put. Long story short, better create a new patch entry instead of using the same one to avoid this kind ofconfusion.

[dmorissette]

Ha! Thanks for the explanation. Trac's behavior with attachments is kind of confusing. I was aware that Trac didn't issue an email when attachments are submitted or modified, but I didn't realize that when an attachment is updated, it just moves the comments in the history as if they were fresh, causing potential confusion.