Monday, March 17, 2008

How to find and use Bursa Wolf parameters for accurate coordinate transformation

Many times I had and heard about the problem of missing Bursa Wolf parameters for more accurate re-projections of data. While GRASS deals with them properly, most of the other GIS and also postgis do not. In Italy we have one of these odd problems related to the fact that historically we used to use as the main projection the EPSG:3003, Monte Mario / Italy Zone 1 projection.
Lately, everyone is migrating for standardization to the EPSG:32632, WGS 84 / UTM zone 32N.

What is the problem in all this?
The fact that the epsg 3004 comes like the following:


PROJCS["Monte Mario / Italy zone 1",
GEOGCS["Monte Mario",
DATUM["Monte_Mario",
SPHEROID["International 1924", 6378388.0, 297.0, AUTHORITY["EPSG","7022"]],
AUTHORITY["EPSG","6265"]],
PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]],
UNIT["degree", 0.017453292519943295],
AXIS["Lon", EAST],
AXIS["Lat", NORTH],
AUTHORITY["EPSG","4265"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["central_meridian", 9.0],
PARAMETER["latitude_of_origin", 0.0],
PARAMETER["scale_factor", 0.9996],
PARAMETER["false_easting", 1500000.0],
PARAMETER["false_northing", 0.0],
UNIT["m", 1.0],
AXIS["x", EAST],
AXIS["y", NORTH],
AUTHORITY["EPSG","3003"]]

which doesn't contain the Bursa Wolf parameters (to know more about it read here) without which there can be reprojection problems of many meters.

So we need a way to get those Bursa Wolf parameters.
Thank God there is a nice site called: http://www.epsg-registry.org

Jump in there and click on the upper tab called: retrieve by code.
You could also search in the query by filter tab, this will be explained later in the post.

I will insert epsg 1660 and click on retrieve, which supplies me the following:



Select your area and click on Report selected results.



Please note that the domain of validity is described. Italy mainland, which is what I need.

Fill in what is asked and just press play :)


On the lower part you will find the needed Bursa Wolf parameters, often also referred to as the "towgs parameters" :)

What is needed to be done, is to add to the 3003 definition these parameters:

TOWGS84[-104.1, -49.1, -9.9, 0.971, -2.917, 0.714, -11.68],

which will result in having the following definition:


PROJCS["Monte Mario / Italy zone 1",
GEOGCS["Monte Mario",
DATUM["Monte Mario",
SPHEROID["International 1924", 6378388.0, 297.0, AUTHORITY["EPSG","7022"]],
TOWGS84[-104.1, -49.1, -9.9, 0.971, -2.917, 0.714, -11.68],
AUTHORITY["EPSG","6265"]],
PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]],
UNIT["degree", 0.017453292519943295],
AXIS["Geodetic longitude", EAST],
AXIS["Geodetic latitude", NORTH],
AUTHORITY["EPSG","4265"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["central_meridian", 9.0],
PARAMETER["latitude_of_origin", 0.0],
PARAMETER["scale_factor", 0.9996],
PARAMETER["false_easting", 1500000.0],
PARAMETER["false_northing", 0.0],
UNIT["m", 1.0],
AXIS["Easting", EAST],
AXIS["Northing", NORTH],
AUTHORITY["EPSG","10003003"]]

Note that I also changed the 3003 to 10003003 to be sure tyo not have conflicts with other crs. I created that in postgis as addionional CRS in order to be able to do proper transformations.

Usually you will not have the proper epsg code to use, right?
Luckily Antonio (thanks man, got a link to share?) supplied us with a better method to find what you are looking for.
Assuming you are still on the epsg-registry site:
1. start by clicking on "query by filter" button;
2. click on "Type" textbox and choose Coordinate Transformation and then Single only;
3. write the name of the area of use (e.g. Italy) in the "Area" textbox and finally click on Search button.
You will get a whole list of possible transformations as a result, which will probably be too many, so if you want to limit them to your area of interest to get the most precisely defined transformation for your region, you have to follow step 1. and 2. the same as before before.
Then:
2. write the name of the transformation of use (e.g. "Monte Mario to WGS 84") in the "Name" textbox
3. click on "Type" textbox and choose Coordinate Transformation, Single only and finally click on Search button.


In JGrass / uDig this is way easier. Click on the button reporting the projection definition. A window will ask to choose the crs from a list. A second tab says: Custom CRS. Select that one and paste your crs definition (in my case the 10003003) into the box and press apply and ok.

Thanks to the GFOSS.it mailinglist for discussing this and to Antonio Falciano for reporting the epsg-registry site.

8 comments:

Anonymous said...

Complimenti e grazie.

Antonio Falciano said...

A simple way to search all transformation parameters on EPSG Geodetic Parameter Dataset

Follow these simple steps:
1. go to http://www.epsg-registry.org/;
2. click on "query by filter" button;
3. click on "Type" textbox and choose Coordinate Transformation and then Single only;
4. write the name of the area of use (e.g. Italy) in the "Area" textbox and finally click on Search button.
The possible transformations here available are too much, so if you want to limit your area of use to a well defined transformation, then you have to:
1. & 2. as showed before;
3. write the name of the transformation of use (e.g. "Monte Mario to WGS 84") in the "Name" textbox;
3. click on "Type" textbox and choose Coordinate Transformation, Single only and finally click on Search button.

Have good transformations!

moovida said...

Antonio, thanks for your great contribution!!
I'm going to update this.

Anonymous said...

Hi Andrea,
I have tried to apply your fantastic instructions to a layer of mine inside jgrass.
I set the TOWGS84 parameter for this layer, but if I open again the CRS window I can not see more TOWGS84 parameter. You can see it here:
http://screencast.com/t/LscNER58

Is it normal?

Thank you,

a

moovida said...

Hi Andrea, I'm not sure about whatyou need to gain, so I try to give you a couple of answers:
1) you are happy with the prj files of your shapes and you want to put them on the map that has the right bursa. In that case you should, before you add layer to the map, set the map's projection to the proper projection and only after that add the maps, so they get reprojected on the fly. If the layer are already there, they are reprojected from their actual situation alltogether the same way.
2) You want to make your shapefiles know the proper (bursa added) projection. In that case the best to do is to add the customized projection WKT directly into the prj file of the shapefile. This method is the less errorprone, since the map is aware of the proper projection.

I hope this was of help. If you have problems on this, feel free to ask. Probably at that piont we could take the discussion to the mailinglist.


PS: nice way to do screencasts. I will try that out and create JGrass tutorials in future :) A pitty that they do not support voice (I guess)

Anonymous said...

I will reply you in the mailing list.

With Jing you can also record voice.

ciao,

Andrea

Anonymous said...

great! it worked very well for me!

Anonymous said...

Sei stato un'aiuto meraviglioso!!
Se potessi ti farei santo oggi!!

Grazie

Carlo