Problem with PLY and REF Points when Georeferencing BSB Charts

[draagedo]

I have an issue when georeferencing png-files, and hope someone here can help me.

I am automatically downloading charts from statkart.no, and try to georeference them and convert them into bsb-format on the fly. The header file (used with bsblib) made for this particular chart is below, and the problem can be seen in the screenshot attached.

VER/2.0
BSB/NA=Et navn
    NU=,RA=4048,2825,DU=300
KNP/SC=50000,GD=,PR=MERCATOR,PP=UNKNOWN
    PI=0.0,SP=Unknown,SK=0.0
    UN=METERS,SD=EQUINOCTIAL SPRING LOW WATER,DX=0.136,DY=0.136
REF/1,1,1,59.100000,5.400000
REF/2,1,2825,58.900000,5.400000
REF/3,4048,1,59.100000,5.900000
REF/4,4048,2825,58.900000,5.900000
PLY/1,58.900000,5.900000
PLY/2,59.100000,5.900000
PLY/3,59.100000,5.400000
PLY/4,58.900000,5.400000
OST/1
DTM/0.0,0.0

As you can see the PLY points and the REF points contain the same data, but the actual chart picture does not fill up the PLY area.

Does anyone have an idea why I get this problem? I suspect the header file myself. The position seams to be quite right with regards to the latitude, but not, as you can guess, the longitude.

[cagney]

Originally Posted by draagedo

I have an issue when georeferencing png-files, and hope someone here can help me.

I am automatically downloading charts from statkart.no, and try to georeference them and convert them into bsb-format on the fly. The header file (used with bsblib) made for this particular chart is below, and the problem can be seen in the screenshot attached.

VER/2.0
BSB/NA=Et navn
    NU=,RA=4048,2825,DU=300
KNP/SC=50000,GD=,PR=MERCATOR,PP=UNKNOWN
    PI=0.0,SP=Unknown,SK=0.0
    UN=METERS,SD=EQUINOCTIAL SPRING LOW WATER,DX=0.136,DY=0.136
REF/1,1,1,59.100000,5.400000
REF/2,1,2825,58.900000,5.400000
REF/3,4048,1,59.100000,5.900000
REF/4,4048,2825,58.900000,5.900000
PLY/1,58.900000,5.900000
PLY/2,59.100000,5.900000
PLY/3,59.100000,5.400000
PLY/4,58.900000,5.400000
OST/1
DTM/0.0,0.0

As you can see the PLY points and the REF points contain the same data, but the actual chart picture does not fill up the PLY area.

Does anyone have an idea why I get this problem? I suspect the header file myself. The position seams to be quite right with regards to the latitude, but not, as you can guess, the longitude.

If your REF points are correct, this probably indicates that your picture is not pure Mercator. To save the situation use more REF points. Try a 4x4 pattern for example, to start with, and check the result.

Thomas

[draagedo]

But all the REF points are in the corners. Should not the corners be "right" then?

About the Meracator you are probably right. It is an UTM projection, but changing the header file does not help. I will have to try adding more REF points then, but I can't see that helping since the corners are off by quite a bit.

[Totobeloeil]

draagedo,

On a windows platform, you can use mapcal_2 to calibrate charts.

Then, KAP/BSB files can be generated by using the following programs :
mc2bsbh, nconvert, tif2bsb.

mapcal_2 comes with seaclear, tif2bsb with libbsb and nconvert with xnview.
Could be easily found on the web.

mc2bsbh was made by dacust :
MapCal Utilities

A cmd tool (kapgen) was made by GPS-Marco :
MCE66 - Software

Take care : this is beta software, made for special purpose.
it is downloading charts before proceeding them.

The following thread could also help you :
http://www.cruisersforum.com/forums/f134/charts-31346-38.html#post404290

Jean-Pascal

[HamishB]

Originally Posted by draagedo

It is an UTM projection

is the original image in a georeferenced format? i.e. does "gdalinfo" know where it is before you calibrate it?

if so, you can use "gdalwarp" to reproject it into plain Mercator from UTM.


Hamish

ps, check that all your x's and y's are not reversed somewhere.

[HamishB]

.. and if they aren't georeferenced you can use gdal_translate to convert the PNG + known bounds & utm zone into a GeoTiff, and then use gdalwarp to reproject it into a custom Mercator with a lat of true scale in the middle of your map somewhere. gdalinfo will then tell you the lat/lon of the four corners.

[draagedo]

The image is not georeferenced. I wrote a script that takes two sets of wgs84 coordinates in, converts them to UTM and downloads the png file

i.e.lower_lat = 58.9
lower_long = 5.4
upper_lat = 59.1
upper_long = 5.9

gives you this link: http://openwms.statkart.no/skwms1/wm...NSPARENT=FALSE

based on the conversion found from pygps and here

Have tried to use chartcal II from the seaclear program, but I don't have a windows machine (running it under wine) so I have not yet tried the conversion tool to a bsb header, but done it manually.

I will give chartcal a few more shots, and georeference from other known locations on the map instead of the corners I have used until now. I still find it strange that the corners from the PLY points are so far from the corners of the calibration REF points when they are the very same coordinates in the header file.

[HamishB]

Originally Posted by draagedo

The image is not georeferenced. I wrote a script that takes two sets of wgs84 coordinates in, converts them to UTM and downloads the png file

i.e.lower_lat = 58.9
lower_long = 5.4
upper_lat = 59.1
upper_long = 5.9

gives you this link: http://openwms.statkart.no/skwms1/wms.sjo_hovedkart2?VERSION=1.1.1&SERVICE=WMS&REQUE ST=GetMap&SRS=EPSG:32632&BBOX=292635.96,6534498.05 ,322455.46,6555310.02&LAYERS=Sjokart_Hovedkartseri en2&STYLES=&WIDTH=4048&HEIGHT=2825&FORMAT=image/png&BGCOLOR=0xFFFFFF&TRANSPARENT=FALSE

Ok, so it downloads it from a WMS server.

with &SRS=EPSG:32632 you control the output map projection. /usr/share/proj/epsg tells me that 32632 is UTM32. (unsurprising )

you can do a capabilities request to see what other projections are available if needed. (XML)

also you might see if GeoTiff is a supported by the WMS server. If you can do
&FORMAT=image/geotiff to get a georeferenced image.

I am surprised it lets you download something this big in a single tile:
&WIDTH=4048&HEIGHT=2825

based on the conversion found from pygps and here

cs2cs from PROJ.4 does the trick for me:

echo "long lat" | cs2cs +init=epsg:4326 +to +init=epsg:32632

I will give chartcal a few more shots, and georeference from other known locations on the map instead of the corners I have used until now. I still find it strange that the corners from the PLY points are so far from the corners of the calibration REF points when they are the very same coordinates in the header file.

The image may be constrained to the realities of DX,DY and image width/height while the PLY red bounding box is not.


Note you can take your known projection (EPSG:32632) and write a little world file to make the PNG georeferenced in the eyes of gdalwarp:
see http://www.gdal.org/frmt_various.html#PNG


Hamish

[brak]

Where is the original image file? If you want to share the georeferenced chart - I'd be interested in taking a look at it.

Originally Posted by draagedo

The image is not georeferenced. I wrote a script that takes two sets of wgs84 coordinates in, converts them to UTM and downloads the png file

i.e.lower_lat = 58.9
lower_long = 5.4
upper_lat = 59.1
upper_long = 5.9

gives you this link: http://openwms.statkart.no/skwms1/wm...NSPARENT=FALSE

based on the conversion found from pygps and here

Have tried to use chartcal II from the seaclear program, but I don't have a windows machine (running it under wine) so I have not yet tried the conversion tool to a bsb header, but done it manually.

I will give chartcal a few more shots, and georeference from other known locations on the map instead of the corners I have used until now. I still find it strange that the corners from the PLY points are so far from the corners of the calibration REF points when they are the very same coordinates in the header file.

[HamishB]

Originally Posted by brak

Where is the original image file?

You pass the request to the WMS server as a HTTP GET, and it responds by sending you the appropriately prepared image.

The so original image file is the URL.


Hamish

[brak]

Ha, indeed it is Thanks, downloading now.

Originally Posted by HamishB

You pass the request to the WMS server as a HTTP GET, and it responds by sending you the appropriately prepared image.

The so original image file is the URL.


Hamish

[brak]

So I used the same KAP header that was provided in the first message with the image, gluing them together using tif2bsb.
The result is a reasonable chart. It works in PolarView and a few other viewers I tested. One issue may be that the chart image is not in Mercator format. That, among other things, results in it having to be scaled differently along vertical and horizontal axis. Specifically, the image has to be "stretched" vertically. Depending on how OpenCPN deals with that, I would suspect that this is the issue (i.e. it looks like the image area is longer, as it should be - but the actual pixels are not stretched to fill it).

Here is what the image looks like when it is stretched to the coverage area. Note the shape of navigation light circles - they are ellipses, stretched vertically. Note that it even matches general outlines of an underlying vector chart, though that chart is of much smaller scale so naturally there is some difference.

Added: yep, that's the issue. The original screen capture shows that light sectors are perfectly round. Looks like it's trying to scale the image equally in all directions. That'll do it.

[draagedo]

Originally Posted by brak

One issue may be that the chart image is not in Mercator format.

That is quite right. I have changed it to Transverse Mercator, but is does not help. I have changed the dx and dy as well as they were wrong.

Originally Posted by brak

That, among other things, results in it having to be scaled differently along vertical and horizontal axis. Specifically, the image has to be "stretched" vertically. Depending on how OpenCPN deals with that, I would suspect that this is the issue (i.e. it looks like the image area is longer, as it should be - but the actual pixels are not stretched to fill it).

I guess I can try and stretch the image myself. I just found openCPN a few weeks ago, and have learned a lot about chart projection since then :P Maybe someone else can fill in how and if OCPN do any stretching?

Originally Posted by brak

Added: yep, that's the issue. The original screen capture shows that light sectors are perfectly round. Looks like it's trying to scale the image equally in all directions. That'll do it.

Can I remedy this in any way?

It looks like the actual projection is not that far off now, and the error I see is in my own location is because of the concrete buliding I am living in. Will have to confirm this with a few positions from a "real" chart...

[cagney]

Originally Posted by draagedo

That is quite right. I have changed it to Transverse Mercator, but is does not help. I have changed the dx and dy as well as they were wrong.


I guess I can try and stretch the image myself. I just found openCPN a few weeks ago, and have learned a lot about chart projection since then :P Maybe someone else can fill in how and if OCPN do any stretching?



Can I remedy this in any way?

It looks like the actual projection is not that far off now, and the error I see is in my own location is because of the concrete buliding I am living in. Will have to confirm this with a few positions from a "real" chart...

OpenCPN only handles Mercator projections, no matter what PR value you give. Of course,native handling of tmerc charts is on the wish list.
Once again, the way to make this chart show up reasonable in OpenCPN, not counting converting from TMerc to Merc with gdalwarp, is to throw in many more REF:s. This forces OpenCPN to use a polynomial solution instead of a "pure" merc solution. The result depends on how skewed the original chart is, in other words how far is this chart from the central meridian used in the tmerc projection. Don't expect a perfect fit, but a "reasonable" one.
Thomas

[HamishB]

Originally Posted by draagedo

I have changed it to Transverse Mercator, but is does not help.

Cagney touches on this, but I think it bears a lot of explaining:

Mercator and Transverse Mercator are not the same thing.

Mercator is like the standard map of the world you see with Greenland distorted to be bigger than Australia. Only areas along the equator (the central parallel) remain undistorted. More local charts using the Mercator projection will reset the standard parallel to be a nearby latitude, but the distortion remains as you move away from it.

Mercator is all but dead as a map projection outside of the nautical/navigation world. (or at least it should be)


Transverse Mercator (which includes the UTMs) is the same but different. Trans. Merc. turns the projection by 90 degrees and uses a central meridian instead of a central parallel (as the word "transverse" implies). So the line north-south along the "lon_0" longitude is undistorted and it gets more and more distorted the further east-west you move from it. If your chart covers more than a couple hundred miles east-west, you really shouldn't be using tmerc.


So applying mercator transform maths to images projected with transverse mercator isn't going to show you what you want, the hourglass shaped correction will be applied to the y axis instead of the x axis.


Hamish