Re: [mkgmap-dev] HGT - getElevation()
Hi,
I'm attaching a next patch for bicubic interpolation. It include
following changes:
- some errors corrected,
- more statistic of interpolation,
- changed condition for applying bicubic interpolation.
In this version bicubic interpolation is applied, when required DEM
resolution is not lower than 1/3 of HGT resolution. I assume, that for
low DEM resolution (big dem-dists values), there is no reason to make
precise interpolation.
I have uploaded compiled mkgmap:
http://files.mkgmap.org.uk/download/404/mkgmap-dem-tdb-r4071-bicubic-6.zip
--
Best regards,
Andrzej
Index: src/uk/me/parabola/imgfmt/app/dem/DEMFile.java
===
--- src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (revision 4071)
+++ src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (working copy)
@@ -52,20 +52,56 @@
* @param outsidePolygonHeight the height value that should be used for
points outside of the bounding polygon
*/
public void calc(Area area, java.awt.geom.Area demPolygonMapUnits,
String pathToHGT, List pointDistances, short outsidePolygonHeight) {
- HGTConverter hgtConverter = new HGTConverter(pathToHGT, area,
demPolygonMapUnits);
+ // HGT area is extended by 0.01 degree in each direction
+ HGTConverter hgtConverter = new HGTConverter(pathToHGT, area,
demPolygonMapUnits, 0.01D);
hgtConverter.setOutsidePolygonHeight(outsidePolygonHeight);
+ int top = area.getMaxLat() * 256;
+ int bottom = area.getMinLat() * 256;
+ int left = area.getMinLong() * 256;
+ int right = area.getMaxLong() * 256;
+
int zoom = 0;
int lastDist = pointDistances.get(pointDistances.size()-1);
for (int pointDist : pointDistances) {
- DEMSection section = new DEMSection(zoom++, area,
hgtConverter, pointDist, pointDist == lastDist);
+ int distance = pointDist;
+ if (distance == -1) {
+ int res = (hgtConverter.getHighestRes() > 0) ?
hgtConverter.getHighestRes() : 1200;
+ distance = (int) Math.round((1 << 29) / (res *
45.0D));
+ }
+ // last 4 bits of distance should be 0
+ distance = ((distance + 8)/16)*16;
+
+ int xtop = top;
+ int xleft = left;
+
+ // align DEM to distance raster, if distance not bigger
than widening of HGT area
+ if (distance < (int)Math.floor((0.01D/45.0D * (1 <<
29 {
+ if (xtop >= 0) {
+ xtop -= xtop % distance;
+ if (xtop < 0x3FFF - distance)
+ xtop += distance;
+ }
+ else {
+ xtop -= xtop % distance;
+ }
+
+ if (xleft >= 0) {
+ xleft -= xleft % distance;
+ }
+ else {
+ xleft -= xleft % distance;
+ if (xleft > Integer.MIN_VALUE +
distance)
+ xleft -= distance;
+ }
+ }
+
+ DEMSection section = new DEMSection(zoom++, xtop,
xleft, xtop - bottom, right - xleft, hgtConverter, distance, pointDist ==
lastDist);
demHeader.addSection(section);
}
return;
}
-
-
public void write() {
ImgFileWriter w = getWriter();
if (w instanceof BufferedImgFileWriter) {
Index: src/uk/me/parabola/imgfmt/app/dem/DEMSection.java
===
--- src/uk/me/parabola/imgfmt/app/dem/DEMSection.java (revision 4071)
+++ src/uk/me/parabola/imgfmt/app/dem/DEMSection.java (working copy)
@@ -48,24 +48,23 @@
private int maxHeight = Integer.MIN_VALUE;
private List tiles = new ArrayList<>();
- public DEMSection(int zoomLevel, Area bbox, HGTConverter hgtConverter,
int pointDist, boolean lastLevel) {
+ public DEMSection(int zoomLevel, int areaTop, int areaLeft, int
areaHeight, int areaWidth, HGTConverter hgtConverter, int pointDist, boolean
lastLevel) {
this.zoomLevel = zoomLevel;
this.lastLevel = lastLevel;
- if (pointDist == -1) {
- int res = (hgtConverter.getHighestRes() > 0) ?
hgtConverter.getHighestRes() : 1200;
-
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, > Bicubic expects grid [-1,0,1,2] X [-1,0,1,2] > fillArray uses [0,1,2,3]X[0,1,2,3] That's the reason for "-1" in following code: h = rdr.ele(xLeft + x - 1, yBottom + y - 1); eleArray[x][y] = h; I have corrected some errors, will post a new patch soon. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi, here is my version of mkgmap for tests: http://files.mkgmap.org.uk/download/403/mkgmap-dem-tdb-r4071-bicubic.zip Changes are following: - dem-dists values are rounded to 16, - DEM coordinates are multiples of dem-dists, - DEM values are calculated with bicubic interpolation for layers with dem-dists <= 9942. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, > Maybe you can create a 2nd branch so that others can try the > different results? I hardly know java, have no IDE for work or any experience in SVN. Basically I can tweak your code and compile it with ant. I could post compiled version if there is any interest. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, I'll try this on monday, Maybe you can create a 2nd branch so that others can try the different results? BTW: I think that the calculation of qx and qy for the interpolation should be changed. The current code in 4071 adds up rounding errors when data for the upper right area is calculated in a large tile spanning multiple degrees. I am working on that. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Samstag, 20. Januar 2018 00:50:19 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have included bicubic interpolation created by Ken Perlin: https://github.com/mlevans/pretty-heatmaps/blob/master/DensityArrayCreation/src/DensityArrayCreation/Bicubic.java Bicubic interpolation is used only for layers with good resolution (3 arc seconds or better). For other cases bilinear interpolation is used. Patch is attached. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi, quick comparison between r4068 original and r4070 with bicubic interpolation. Both use pbf aligned to HGT and the same 3 arc degree HGT with many voids. http://files.mkgmap.org.uk/download/401/orig-aligned.png http://files.mkgmap.org.uk/download/402/bicu-aligned.png DEM shading moved again. DEM heights are similar on both maps, there is difference like 20-30m between DEM peaks. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
I have included bicubic interpolation created by Ken Perlin:
https://github.com/mlevans/pretty-heatmaps/blob/master/DensityArrayCreation/src/DensityArrayCreation/Bicubic.java
Bicubic interpolation is used only for layers with good resolution (3
arc seconds or better). For other cases bilinear interpolation is used.
Patch is attached.
--
Best regards,
Andrzej
Index: src/uk/me/parabola/imgfmt/app/dem/DEMFile.java
===
--- src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (revision 4070)
+++ src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (working copy)
@@ -52,20 +52,57 @@
* @param outsidePolygonHeight the height value that should be used for
points outside of the bounding polygon
*/
public void calc(Area area, java.awt.geom.Area demPolygonMapUnits,
String pathToHGT, List pointDistances, short outsidePolygonHeight) {
- HGTConverter hgtConverter = new HGTConverter(pathToHGT, area,
demPolygonMapUnits);
+ // HGT area is extended by 0.01 degree in each direction
+ HGTConverter hgtConverter = new HGTConverter(pathToHGT, area,
demPolygonMapUnits, 0.01D);
hgtConverter.setOutsidePolygonHeight(outsidePolygonHeight);
+ int top = area.getMaxLat() * 256;
+ int bottom = area.getMinLat() * 256;
+ int left = area.getMinLong() * 256;
+ int right = area.getMaxLong() * 256;
+
int zoom = 0;
int lastDist = pointDistances.get(pointDistances.size()-1);
for (int pointDist : pointDistances) {
- DEMSection section = new DEMSection(zoom++, area,
hgtConverter, pointDist, pointDist == lastDist);
+ int distance = pointDist;
+ if (distance == -1) {
+ int res = (hgtConverter.getHighestRes() > 0) ?
hgtConverter.getHighestRes() : 1200;
+ distance = (int) ((1 << 29) / (res * 45));
+ }
+ // last 4 bits of distance should be 0
+ distance = ((distance + 8)/16)*16;
+
+ int xtop = top;
+ int xleft = left;
+
+ // align DEM to distance raster, if distance not bigger
than widening of HGT area
+ if (distance < (int)Math.floor((0.01D/45.0D * (1 <<
29 {
+ if (xtop >= 0) {
+ xtop -= xtop % distance;
+ if (xtop < 0x3FFF - distance)
+ xtop += distance;
+ }
+ else {
+ xtop -= xtop % distance;
+ }
+
+ if (xleft >= 0) {
+ xleft -= xleft % distance;
+ }
+ else {
+ xleft -= xleft % distance;
+ if (xleft > Integer.MIN_VALUE +
distance)
+ xleft -= distance;
+ }
+ }
+
+ // to improve: DEM secition should get area in 32bit
units
+ DEMSection section = new DEMSection(zoom++, xtop,
xleft, xtop - bottom, right - xleft, hgtConverter, distance, pointDist ==
lastDist);
demHeader.addSection(section);
}
return;
}
-
-
public void write() {
ImgFileWriter w = getWriter();
if (w instanceof BufferedImgFileWriter) {
Index: src/uk/me/parabola/imgfmt/app/dem/DEMSection.java
===
--- src/uk/me/parabola/imgfmt/app/dem/DEMSection.java (revision 4070)
+++ src/uk/me/parabola/imgfmt/app/dem/DEMSection.java (working copy)
@@ -48,16 +48,12 @@
private int maxHeight = Integer.MIN_VALUE;
private List tiles = new ArrayList<>();
- public DEMSection(int zoomLevel, Area bbox, HGTConverter hgtConverter,
int pointDist, boolean lastLevel) {
+ public DEMSection(int zoomLevel, int areaTop, int areaLeft, int
areaHeight, int areaWidth, HGTConverter hgtConverter, int pointDist, boolean
lastLevel) {
this.zoomLevel = zoomLevel;
this.lastLevel = lastLevel;
- if (pointDist == -1) {
- int res = (hgtConverter.getHighestRes() > 0) ?
hgtConverter.getHighestRes() : 1200;
- pointDist = (int) ((1 << 29) / (res * 45));
- }
- this.top = bbox.getMaxLat() * 256;
-
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, > I'd prefer to use unmodified dem-dist values for now, besides that the > unit tests don't compile with it. I'm trying to get DEM looking like in Garmins img. The actual rounding in code can be disabled, but then it would be up to user, to set correct distances. I don't know, if rounding is important. Maybe some devices need it? > After adding 0.01 as new parm some tests fail. This looked like an easy way to get HGT without computing and comparing values for all layers. It could be done better, if the idea fo this patch worked out. Sorry about tests, I don't know how to execute them. Now I consider to add bicubic interpolation, like this code: https://github.com/mlevans/pretty-heatmaps/blob/master/DensityArrayCreation/src/DensityArrayCreation/Bicubic.java -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, I'd prefer to use unmodified dem-dist values for now, besides that the unit tests don't compile with it. After adding 0.01 as new parm some tests fail. Probably not so important, as the unit test HGTConverterTest is bad anyway, it tests non-public methods and your patch shows that this is a bad idea. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Freitag, 19. Januar 2018 02:59:14 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have prepared patch according to these guessed rules. It rounds dem-dist and shift coordinates of top-left corner of DEM areas. Since dem-dist doesn't correspond exactly to HGT, interpolation is always active. I think there could be 2 ways to improve DEM quality. Interpolation could be changed to bicubic or mkgmap could use preprocessed data with correct pixel size. In the latter case probably data format other than HGT should be used. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
I have prepared patch according to these guessed rules. It rounds
dem-dist and shift coordinates of top-left corner of DEM areas. Since
dem-dist doesn't correspond exactly to HGT, interpolation is always active.
I think there could be 2 ways to improve DEM quality. Interpolation
could be changed to bicubic or mkgmap could use preprocessed data with
correct pixel size. In the latter case probably data format other than
HGT should be used.
--
Best regards,
Andrzej
Index: src/uk/me/parabola/imgfmt/app/dem/DEMFile.java
===
--- src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (revision 4070)
+++ src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (working copy)
@@ -52,20 +52,57 @@
* @param outsidePolygonHeight the height value that should be used for
points outside of the bounding polygon
*/
public void calc(Area area, java.awt.geom.Area demPolygonMapUnits,
String pathToHGT, List pointDistances, short outsidePolygonHeight) {
- HGTConverter hgtConverter = new HGTConverter(pathToHGT, area,
demPolygonMapUnits);
+ // HGT area is extended by 0.01 degree in each direction
+ HGTConverter hgtConverter = new HGTConverter(pathToHGT, area,
demPolygonMapUnits, 0.01D);
hgtConverter.setOutsidePolygonHeight(outsidePolygonHeight);
+ int top = area.getMaxLat() * 256;
+ int bottom = area.getMinLat() * 256;
+ int left = area.getMinLong() * 256;
+ int right = area.getMaxLong() * 256;
+
int zoom = 0;
int lastDist = pointDistances.get(pointDistances.size()-1);
for (int pointDist : pointDistances) {
- DEMSection section = new DEMSection(zoom++, area,
hgtConverter, pointDist, pointDist == lastDist);
+ int distance = pointDist;
+ if (distance == -1) {
+ int res = (hgtConverter.getHighestRes() > 0) ?
hgtConverter.getHighestRes() : 1200;
+ distance = (int) ((1 << 29) / (res * 45));
+ }
+ // last 4 bits of distance should be 0
+ distance = ((distance + 8)/16)*16;
+
+ int xtop = top;
+ int xleft = left;
+
+ // align DEM to distance raster, if distance not bigger
than widening of HGT area
+ if (distance < (int)Math.floor((0.01D/45.0D * (1 <<
29 {
+ if (xtop >= 0) {
+ xtop -= xtop % distance;
+ if (xtop < 0x3FFF - distance)
+ xtop += distance;
+ }
+ else {
+ xtop -= xtop % distance;
+ }
+
+ if (xleft >= 0) {
+ xleft -= xleft % distance;
+ }
+ else {
+ xleft -= xleft % distance;
+ if (xleft > Integer.MIN_VALUE +
distance)
+ xleft -= distance;
+ }
+ }
+
+ // to improve: DEM secition should get area in 32bit
units
+ DEMSection section = new DEMSection(zoom++, xtop,
xleft, xtop - bottom, right - xleft, hgtConverter, distance, pointDist ==
lastDist);
demHeader.addSection(section);
}
return;
}
-
-
public void write() {
ImgFileWriter w = getWriter();
if (w instanceof BufferedImgFileWriter) {
Index: src/uk/me/parabola/imgfmt/app/dem/DEMSection.java
===
--- src/uk/me/parabola/imgfmt/app/dem/DEMSection.java (revision 4070)
+++ src/uk/me/parabola/imgfmt/app/dem/DEMSection.java (working copy)
@@ -48,16 +48,12 @@
private int maxHeight = Integer.MIN_VALUE;
private List tiles = new ArrayList<>();
- public DEMSection(int zoomLevel, Area bbox, HGTConverter hgtConverter,
int pointDist, boolean lastLevel) {
+ public DEMSection(int zoomLevel, int areaTop, int areaLeft, int
areaHeight, int areaWidth, HGTConverter hgtConverter, int pointDist, boolean
lastLevel) {
this.zoomLevel = zoomLevel;
this.lastLevel = lastLevel;
- if (pointDist == -1) {
- int res = (hgtConverter.getHighestRes() > 0) ?
hgtConverter.getHighestRes() : 1200;
- pointDist = (int) ((1
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I have executed your script and looked at maps. I have tested them with artificial HGT too. Only the original-unaligned map has DEM at wrong position. It is moved by about 0.00033 degree right and down. I looked at offset calculated as (DEM coordinate)%(dem-dist) and only for this map, offset was bigger than half of dem-dist. Any attempt to align HGT with data reduce this offset, but it still remain quite big, like 1300-3600, while dem-dist is 9942. My guess is, that Mapsource/BaseCamp calculate DEM pixel positions globally, basing on dem-dist. If coordinates of DEM layers aren't aligned with with dem-dist and misalignment is too big, than Mapsource can fetch wrong DEM values. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, these are good hints. I already wondered why upper left corner is not the same at different zoom levels. I'll try to change the code to match these conditions, maybe this solves the problem with "dem movement". Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Donnerstag, 18. Januar 2018 03:08:48 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have looked a bit at Garmin's DEM and I have some guesses. * Coordinate resolution for DEM is 28 bits. Since it is shifted to the left, 4 lower bits are always 0. * The same goes for pixel distance, lower 4 bits should be 0. This is the reason, why Garmin uses 3312 (CF0) instead of 3314 (CF2). * Coordinates of left top corner for each DEM layer are always exact multiple of pixel distance. Since layers get different distances, each layer get a bit different coordinates too. I would expect, that aligning of HGT in mkgmap should at least support 3-rd requirement, but it is not the case. I don't know why. Differences looks too big to be explained by rounding errors. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, I cannot yet confirm this, please see http://files.mkgmap.org.uk/download/400/demmove.7z It contains a script to create different maps and screenshots that show that maps is changed with your patch Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Mittwoch, 17. Januar 2018 23:24:35 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have tried to verify, if aligned DEM is correct. I think algorithm have passed positively 2 tests: - reading height from map gives correct values, - changing bbox of source data doesn't change resulting map, even if it changes offset applied for aligning. The aligning changes a bit shading, but it is difficult to decide which version better. Only your map of Bolivia suggested, that nonaligned version could be better. Maybe we should investigate this one? For example we can oversample HGT using some good algorithm, like spline or lanczos, make unaligned map with oversampled HGT and compare changes in shadings? -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I have looked a bit at Garmin's DEM and I have some guesses. * Coordinate resolution for DEM is 28 bits. Since it is shifted to the left, 4 lower bits are always 0. * The same goes for pixel distance, lower 4 bits should be 0. This is the reason, why Garmin uses 3312 (CF0) instead of 3314 (CF2). * Coordinates of left top corner for each DEM layer are always exact multiple of pixel distance. Since layers get different distances, each layer get a bit different coordinates too. I would expect, that aligning of HGT in mkgmap should at least support 3-rd requirement, but it is not the case. I don't know why. Differences looks too big to be explained by rounding errors. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I have tried to verify, if aligned DEM is correct. I think algorithm have passed positively 2 tests: - reading height from map gives correct values, - changing bbox of source data doesn't change resulting map, even if it changes offset applied for aligning. The aligning changes a bit shading, but it is difficult to decide which version better. Only your map of Bolivia suggested, that nonaligned version could be better. Maybe we should investigate this one? For example we can oversample HGT using some good algorithm, like spline or lanczos, make unaligned map with oversampled HGT and compare changes in shadings? -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, I only see movement when I compare a non-aligned tile compiled with unpatched r4068 with a non-aligned tile compiled with a patched one. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Mittwoch, 17. Januar 2018 22:36:08 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have compiled my 2 tiles and the one you found as a problematic. Now with SRTM 3 seconds. I have used mkgmap 4068 with your patch. It all looks good. Shading is the same on all maps, peak DEM height for a mountain is the same and nearly on the same position, which is at HGT node. I don't compare to non-aligned/interpolated versions. I have made screenshots, I can upload them, if you are interested. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I have compiled my 2 tiles and the one you found as a problematic. Now with SRTM 3 seconds. I have used mkgmap 4068 with your patch. It all looks good. Shading is the same on all maps, peak DEM height for a mountain is the same and nearly on the same position, which is at HGT node. I don't compare to non-aligned/interpolated versions. I have made screenshots, I can upload them, if you are interested. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I have prepared artificial HGT with both resolution and tested both. I haven't tested with real HGT 3". I will check this too. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, note that I used 3'' hgt files. Did you try that as well? Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Mittwoch, 17. Januar 2018 20:21:35 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, > So I tried other bboxes, and found this one: > 2902: 2291022,918105 to 2297489,946049 I have compiled this tile and found no difference, it looks the same as my tiles. DEM heights are the same. I have prepared dummy HGT with checkerboard patter. It makes quite easy to verify, if DEM is aligned correctly. DEM contains heights 10m and 20m on a squares with size 0.01*0.01 degree. http://files.mkgmap.org.uk/download/399/test_hgt.7z -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, > So I tried other bboxes, and found this one: > 2902: 2291022,918105 to 2297489,946049 I have compiled this tile and found no difference, it looks the same as my tiles. DEM heights are the same. I have prepared dummy HGT with checkerboard patter. It makes quite easy to verify, if DEM is aligned correctly. DEM contains heights 10m and 20m on a squares with size 0.01*0.01 degree. http://files.mkgmap.org.uk/download/399/test_hgt.7z -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, I've tried your areas.list with 3'' hgt files with three different mkgmap versions (r4068, r4068 + dem-align2.patch (4068_e) and r4068+dem-move2.patch (4068_m) . I looked at 49.20 20.0 which - I think - should give 1872 m. r4068 with aligned pbf shows 1872 m r4068 with unaligned pbf shows 1863 m (due to interpolation) r4068e with aligned pbf shows 1872 m r4068e with unaligned pbf shows 1872m r4068m with aligned pbf shows 1872 m r4068m with unaligned pbf shows 1872m I also did not see any movement in DEM (good), so I wondered why results differ so much from my tests with a tile in bolivia. So I tried other bboxes, and found this one: 2902: 2291022,918105 to 2297489,946049 When you use this for splitter it produces an OSM file with This seems to be a worser case as it shows different results: r4068 says height is 1862, the others say 1872 but this time I see a clear movement in DEM between unpatched and patched binary. Also this one shows this problem: 2900: 2291022,918106 to 2297527,946049 It seems that there is a range where aligning is not okay. Still trying to find out the threshold value(s). Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Dienstag, 16. Januar 2018 22:30:21 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have done 2 kind of tests, both with mkgmap patched with your patch for aligned HGT. Fist test is more or less the same what you propose. I looked for a round peak on a map and then searched for the highest value of DEM on a peak, moving cursor around. With zoom 20m, I tried to estimate area, where DEM is the highest and then I noted coordinates of the center of this area. Then I loaded original HGT into QGIS and looked for the coordinates and position inside HGT pixel, which is displayed by QGIS as a square. I got very good match, coordinates were near the center of the pixel and height was the same as in Mapsource. In second test I calculated manually areas.list with 2 tiles for splitter. One tile aligned with HGT raster, second with left and top edge moved by about half of HGT pixel size. I had to tweak these tiles a bit and the final result was this: 2900: 2291022,918087 to 2297546,946049 # : aligned 2901: 2291022,918097 to 2297534,946049 # : non-aligned After compiling these tiles, I got 2 maps with following coordinates: N: 49.299989, S: 49.15, W: 19.700010, E: 20.39 DEM layers: 1, N: 49.300010, W: 19.700010 N: 49.299731, S: 49.15, W: 19.700224, E: 20.39 DEM layers: 1, N: 49.300010, W: 19.700010 As you can see, DEM position was the same, but map position was a bit different in both tiles. Then I saved both tiles on virtual drive, where they could be seen by BaseCamp as independent maps. I could switch between tiles and compare shading. I uploaded screenshots. I think about 2 more ways of testing. Frank in his manual about DEM, described a clever way to sample DEM values from Mapsource. His tools could be used to validate DEM conversion. Another way could be to artificially move DEM area by multiple of HGT pixel size. This would create a DEM with bigger offset but with the same heights taken from interpolation. If the look of the map wouldn't change, then we could assume, that programs can deal with offset correctly. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I have done 2 kind of tests, both with mkgmap patched with your patch for aligned HGT. Fist test is more or less the same what you propose. I looked for a round peak on a map and then searched for the highest value of DEM on a peak, moving cursor around. With zoom 20m, I tried to estimate area, where DEM is the highest and then I noted coordinates of the center of this area. Then I loaded original HGT into QGIS and looked for the coordinates and position inside HGT pixel, which is displayed by QGIS as a square. I got very good match, coordinates were near the center of the pixel and height was the same as in Mapsource. In second test I calculated manually areas.list with 2 tiles for splitter. One tile aligned with HGT raster, second with left and top edge moved by about half of HGT pixel size. I had to tweak these tiles a bit and the final result was this: 2900: 2291022,918087 to 2297546,946049 # : aligned 2901: 2291022,918097 to 2297534,946049 # : non-aligned After compiling these tiles, I got 2 maps with following coordinates: N: 49.299989, S: 49.15, W: 19.700010, E: 20.39 DEM layers: 1, N: 49.300010, W: 19.700010 N: 49.299731, S: 49.15, W: 19.700224, E: 20.39 DEM layers: 1, N: 49.300010, W: 19.700010 As you can see, DEM position was the same, but map position was a bit different in both tiles. Then I saved both tiles on virtual drive, where they could be seen by BaseCamp as independent maps. I could switch between tiles and compare shading. I uploaded screenshots. I think about 2 more ways of testing. Frank in his manual about DEM, described a clever way to sample DEM values from Mapsource. His tools could be used to validate DEM conversion. Another way could be to artificially move DEM area by multiple of HGT pixel size. This would create a DEM with bigger offset but with the same heights taken from interpolation. If the look of the map wouldn't change, then we could assume, that programs can deal with offset correctly. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, can you explain in detail what you did? Did you use patched mkgmap for this? reg. testing: My current thinking is that we should look at the height values displayed at the status line when in highest zoom (20m). For example, create a map that contains an area around -20.0 -64.0 , use Strg+T to jump to exactly S20 W64 , make sure that the cursor is in the crosshair and it should display 1284m because that is the value that I find as 1st value in s21w064.hgt (it starts with 0x0504). Repeat this procedure with different maps, one with bbox aligned, one with a large difference (1/2 hgt pixel width, maybe also 1/3 and 3/4 or so). All this with unpatched and patched r4068. I assume that we get the best result with an input file that has a bounds statement which places the upper left corner as close as possible to the hgt raster. Does that make sense? Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Dienstag, 16. Januar 2018 16:36:22 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have tried to align pbf bbox with HGT. I have ceated 2 maps, one with god alignment and second with left and top edge moved by 1/2 HGT pixel size. I have used HGT with resolution 3600, so the movement was about 15m. I think shading is the same on both maps. There is no obvious change, but I can't make BaseCamp to display maps at exactly the same position, so there is a small movement of all objects. See pictures: http://files.mkgmap.org.uk/download/396/hgt-aligned.png http://files.mkgmap.org.uk/download/397/hgt-nonaligned.png DEM layer is the same for both maps, position and size exactly the same. I have found, why shading on your map looks more detailed. You have created only single layer of DEM, while I use 4 layers. Shading with multiple layers looks softer. Another weird effect is that with single layer, shading gets softer, when you set "Map Detail" to higher in Mapsource. With multilayer is opposite, probably Mapsource switches between layers. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, sorry, maybe I made again something wrong while testing :-(( This testing is really difficult, maybe I've exchanged results of patched and unpatched version. In this case my previous post was completely wrong. Or one of my unit tests is wrong. I'll check this again tomorrow, need some rest now ;-) Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Gerd Petermann <gpetermann_muenc...@hotmail.com> Gesendet: Dienstag, 16. Januar 2018 17:40:54 An: Development list for mkgmap Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Andrzej, reg. screen shots: I created my two screenshots this way: 1) locate a position -> create 1st screenshot -> leave MapSource 2) replace map content (I just rename dirs, both maps are created with the same options, just different binaries) 3) open MapSource, press Strg+G two times, create 2nd screenshot I think the same works with Basecamp, both remember the position in the map in this case. You can also use Strg+T to jump to a position. I've also tested quality again now and found that DEM data is clearly better when DEM bbox is aligned to hgt raster. So, yes, we should align, and I was wrong, there is no need to change splitter for that, we just have to align the DEM bbox :-) This allows to get the best result out of a hgt file without using a dem-dist value that is lower than the resolution of the hgt file. Question is if anything is improved with e.g. dem-dists=3312,... when dem-dists=9942 gives the best possible result? Gerd P.S. please review new patch, last one did not use moved value for right boundary Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Dienstag, 16. Januar 2018 16:36:22 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have tried to align pbf bbox with HGT. I have ceated 2 maps, one with god alignment and second with left and top edge moved by 1/2 HGT pixel size. I have used HGT with resolution 3600, so the movement was about 15m. I think shading is the same on both maps. There is no obvious change, but I can't make BaseCamp to display maps at exactly the same position, so there is a small movement of all objects. See pictures: http://files.mkgmap.org.uk/download/396/hgt-aligned.png http://files.mkgmap.org.uk/download/397/hgt-nonaligned.png DEM layer is the same for both maps, position and size exactly the same. I have found, why shading on your map looks more detailed. You have created only single layer of DEM, while I use 4 layers. Shading with multiple layers looks softer. Another weird effect is that with single layer, shading gets softer, when you set "Map Detail" to higher in Mapsource. With multilayer is opposite, probably Mapsource switches between layers. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, reg. screen shots: I created my two screenshots this way: 1) locate a position -> create 1st screenshot -> leave MapSource 2) replace map content (I just rename dirs, both maps are created with the same options, just different binaries) 3) open MapSource, press Strg+G two times, create 2nd screenshot I think the same works with Basecamp, both remember the position in the map in this case. You can also use Strg+T to jump to a position. I've also tested quality again now and found that DEM data is clearly better when DEM bbox is aligned to hgt raster. So, yes, we should align, and I was wrong, there is no need to change splitter for that, we just have to align the DEM bbox :-) This allows to get the best result out of a hgt file without using a dem-dist value that is lower than the resolution of the hgt file. Question is if anything is improved with e.g. dem-dists=3312,... when dem-dists=9942 gives the best possible result? Gerd P.S. please review new patch, last one did not use moved value for right boundary Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Dienstag, 16. Januar 2018 16:36:22 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have tried to align pbf bbox with HGT. I have ceated 2 maps, one with god alignment and second with left and top edge moved by 1/2 HGT pixel size. I have used HGT with resolution 3600, so the movement was about 15m. I think shading is the same on both maps. There is no obvious change, but I can't make BaseCamp to display maps at exactly the same position, so there is a small movement of all objects. See pictures: http://files.mkgmap.org.uk/download/396/hgt-aligned.png http://files.mkgmap.org.uk/download/397/hgt-nonaligned.png DEM layer is the same for both maps, position and size exactly the same. I have found, why shading on your map looks more detailed. You have created only single layer of DEM, while I use 4 layers. Shading with multiple layers looks softer. Another weird effect is that with single layer, shading gets softer, when you set "Map Detail" to higher in Mapsource. With multilayer is opposite, probably Mapsource switches between layers. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev dem-move2.patch Description: dem-move2.patch ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I have tried to align pbf bbox with HGT. I have ceated 2 maps, one with god alignment and second with left and top edge moved by 1/2 HGT pixel size. I have used HGT with resolution 3600, so the movement was about 15m. I think shading is the same on both maps. There is no obvious change, but I can't make BaseCamp to display maps at exactly the same position, so there is a small movement of all objects. See pictures: http://files.mkgmap.org.uk/download/396/hgt-aligned.png http://files.mkgmap.org.uk/download/397/hgt-nonaligned.png DEM layer is the same for both maps, position and size exactly the same. I have found, why shading on your map looks more detailed. You have created only single layer of DEM, while I use 4 layers. Shading with multiple layers looks softer. Another weird effect is that with single layer, shading gets softer, when you set "Map Detail" to higher in Mapsource. With multilayer is opposite, probably Mapsource switches between layers. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, okay, got it. You are right, one should test the encoded values. OSM data might as well be shifted due to wrong sat images. I'm still coding the unit tests, next I'll look again at the different versions of the patch. ciao, Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Montag, 15. Januar 2018 23:33:03 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I use SRTM3 v3, so map can be different, but I'm surprised by difference in details. That's why I wanted to recompile map. I can see shift in shading when aligning DEM, I told so in my first post. But I rather don't notice that shading doesn't fit map features. Probably this becomes visible, when comparing 2 versions. Judging shadings can be a bit subjective. I have tried to verify actual height displayed under cursor versus HGT data and this looked good. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I use SRTM3 v3, so map can be different, but I'm surprised by difference in details. That's why I wanted to recompile map. I can see shift in shading when aligning DEM, I told so in my first post. But I rather don't notice that shading doesn't fit map features. Probably this becomes visible, when comparing 2 versions. Judging shadings can be a bit subjective. I have tried to verify actual height displayed under cursor versus HGT data and this looked good. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, hgt files are 3'' from viewfinder, the pictures show an area in Bolivia. The picture shows the location in the status line. And yes, the extended DEM data is clearly shifted relative to the OSM data. I assume now that Garmin uses this feature to match existing map data with existing DEM data. Do you want to say that you don't see a movement in your tests? Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Montag, 15. Januar 2018 21:41:48 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, on your pictures it looks like shading form extended DEM doesn't match map features. Is it HGT with resolution 1200? Shading looks more detailed than on my map, which uses DEM created by BuidDemFile. I have to recompile my map to compare. I have done different test. I have analyzed a round peak on a map and traced coordinates of the center of highest DEM value. All at zoom 20m. Than I looked up the coordinates on HGT tile in QGIS. There were good match, coordinates were max 5m from the center of a HGT pixel. Then I checked offset of DEM, it was about 7m north and 18m west. For me it looks like shifting DEM works correctly, but I can't explain effects on your map. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, on your pictures it looks like shading form extended DEM doesn't match map features. Is it HGT with resolution 1200? Shading looks more detailed than on my map, which uses DEM created by BuidDemFile. I have to recompile my map to compare. I have done different test. I have analyzed a round peak on a map and traced coordinates of the center of highest DEM value. All at zoom 20m. Than I looked up the coordinates on HGT tile in QGIS. There were good match, coordinates were max 5m from the center of a HGT pixel. Then I checked offset of DEM, it was about 7m north and 18m west. For me it looks like shifting DEM works correctly, but I can't explain effects on your map. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, attached patch moves hgt area. With this patch I see two effects: 1) DEM tiles have the same size (values in header), but different data. The aligned tiles are always much bigger, so I think I underestimated the nnegative effect of interpolation. 2) In Mapsource data DEM data is a bit clearer, but it is shifted relative to the other data. This also happens with your patch, no difference visible. See http://files.mkgmap.org.uk/download/394/extend.png and http://files.mkgmap.org.uk/download/393/orig.png If you open both screenshots and switch between them the move is visible. My conclusion: 1) aligning the tile to hgt raster would be better if quality matters 2) if we do that, we must align all data, not only DEM, and I think that means that we should do it in splitter. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Gerd Petermann <gpetermann_muenc...@hotmail.com> Gesendet: Montag, 15. Januar 2018 14:57:17 An: Development list for mkgmap Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Andrzej, ooops , just noticed that I did not update my patch. See new version attached. I am now playing with a version that moves all corners into the same direction. BTW: I've used --dem-poly=bolivia.poly with file from geofabrik. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Montag, 15. Januar 2018 12:51:00 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I guess your code already makes sure, that whole area is covered by DEM. So yes, extending bottom/right is not necessary. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev dem-move.patch Description: dem-move.patch ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, ooops , just noticed that I did not update my patch. See new version attached. I am now playing with a version that moves all corners into the same direction. BTW: I've used --dem-poly=bolivia.poly with file from geofabrik. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Montag, 15. Januar 2018 12:51:00 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I guess your code already makes sure, that whole area is covered by DEM. So yes, extending bottom/right is not necessary. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev dem-align-2c.patch Description: dem-align-2c.patch ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd On 15.01.2018 15:05, Gerd Petermann wrote: > Wouldn't it be better to change splitter so that it aligns tiles to HGT > raster? As long as we have rectangle map tiles, it would be a possible solution. But in future we might have non-rectangle tiles or even tiles, not splitted by splitter. So I think in the end mkgmap should be handle it. Don't know if possible, but maybe mkgmap can virtually extend the tile to the next matching hgt-point. Henning ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I guess your code already makes sure, that whole area is covered by DEM. So yes, extending bottom/right is not necessary. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, I think your patch increases the box too much. I think we only need an aligned upper left corner. See my attached version. Still, even with my patch DEM size can increase quite a lot. I've created a map of bolivia with without --dem-dists option. Size of gmap folder: 223.511k with r4061 238.940k with r4061 + your patch 238.882k with r4061 + my patch So, we have +15 MB for DEM with this alignment. I assume that the increase depends on the areas.list, so maybe my one is a worse case. Wouldn't it be better to change splitter so that it aligns tiles to HGT raster? Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Sonntag, 14. Januar 2018 16:46:29 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, > you just want to achive that the calculated postions are as close as > possible to the positions in the hgt file Yes. Linear interpolation is averaging too. If DEM can use exact values form HGT, then there will be no averaging. Even interpolated values near HGT node should be the same, because of rounding to 1m. I have limited aligning to resolution better or equal 1200 and prepared HGT area for worst case, which would be 1200. Patch attached. -- Best regards, Andrzej dem-align-2b.patch Description: dem-align-2b.patch areas.list Description: areas.list ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
> you just want to achive that the calculated postions are as close as
> possible to the positions in the hgt file
Yes. Linear interpolation is averaging too. If DEM can use exact values
form HGT, then there will be no averaging. Even interpolated values near
HGT node should be the same, because of rounding to 1m.
I have limited aligning to resolution better or equal 1200 and prepared
HGT area for worst case, which would be 1200. Patch attached.
--
Best regards,
Andrzej
Index: src/uk/me/parabola/imgfmt/app/dem/DEMFile.java
===
--- src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (revision 4056)
+++ src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (working copy)
@@ -22,6 +22,7 @@
import uk.me.parabola.imgfmt.app.ImgFileWriter;
import uk.me.parabola.imgfmt.fs.ImgChannel;
import uk.me.parabola.mkgmap.reader.hgt.HGTConverter;
+import uk.me.parabola.imgfmt.Utils;
/**
* The DEM file. This consists of information about elevation. It is used for
hill shading
@@ -52,9 +53,43 @@
* @param outsidePolygonHeight the height value that should be used for
points outside of the bounding polygon
*/
public void calc(Area area, java.awt.geom.Area demPolygonMapUnits,
String pathToHGT, List pointDistances, short outsidePolygonHeight) {
- HGTConverter hgtConverter = new HGTConverter(pathToHGT, area,
demPolygonMapUnits);
+ double top = Utils.toDegrees(area.getMaxLat());
+ double bottom = Utils.toDegrees(area.getMinLat());
+ double left = Utils.toDegrees(area.getMinLong());
+ double right = Utils.toDegrees(area.getMaxLong());
+
+ // expand HGT area for aligning with resolution 1200
+ double hgtDis = 1.0D/1200;
+ double hgtTop = Math.ceil(top/hgtDis)*hgtDis;
+ double hgtBottom = Math.floor(bottom/hgtDis)*hgtDis;
+ double hgtLeft = Math.floor(left/hgtDis)*hgtDis;
+ double hgtRight = Math.ceil(right/hgtDis)*hgtDis;
+ Area hgtArea = new Area(hgtBottom, hgtLeft, hgtTop, hgtRight);
+
+ HGTConverter hgtConverter = new HGTConverter(pathToHGT,
hgtArea, demPolygonMapUnits);
hgtConverter.setOutsidePolygonHeight(outsidePolygonHeight);
+ int res = hgtConverter.getHighestRes();
+ if (res >= 1200) {
+ if (pointDistances.get(0) == -1 ||
pointDistances.get(0) == (int) ((1 << 29) / (res * 45))) {
+ // align area with HGT raster
+ if (res == 1200) {
+ // already calculated
+ area = hgtArea;
+ }
+ else {
+ hgtDis = 1.0D/res;
+ top = Math.ceil(top/hgtDis)*hgtDis;
+ bottom =
Math.floor(bottom/hgtDis)*hgtDis;
+ left = Math.floor(left/hgtDis)*hgtDis;
+ right = Math.ceil(right/hgtDis)*hgtDis;
+
+ area = new Area(bottom, left, top,
right);
+ }
+ System.out.println("DEM aligned to HGT,
resolution " + res);
+ }
+ }
+
int zoom = 0;
int lastDist = pointDistances.get(pointDistances.size()-1);
for (int pointDist : pointDistances) {
___
mkgmap-dev mailing list
mkgmap-dev@lists.mkgmap.org.uk
http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, If I got that right you just want to achive that the calculated postions are as close as possible to the positions in the hgt file. If you request a point that is not within the original area I would expect trouble (either NPE or wrong height 0). If that can happen, you should create a new hgtConverter. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Sonntag, 14. Januar 2018 16:10:36 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, > why you don't use the modified area for the hgtConverter as well? I need HGT resolution first. I guess resolution comes after initialization of hgtConverter. Probably area for hgtConverter should be expanded, or maybe initialized with a bit bigger area? -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, > why you don't use the modified area for the hgtConverter as well? I need HGT resolution first. I guess resolution comes after initialization of hgtConverter. Probably area for hgtConverter should be expanded, or maybe initialized with a bit bigger area? -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, did not try the patch but I wonder why you don't use the modified area for the hgtConverter as well? Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Sonntag, 14. Januar 2018 14:29:57 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have attempted to reduce interpolation, when DEM pixel size is the same as HGT pixel. I have stretched a bit area of DEM, so edges are aligned with HGT. I'm attaching a patch. Resulting map looks similar, only shading has moved a bit. I'm not sure if this is the result of differences in interpolation or maybe change introduced by moving DEM borders. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
I have attempted to reduce interpolation, when DEM pixel size is the
same as HGT pixel. I have stretched a bit area of DEM, so edges are
aligned with HGT. I'm attaching a patch.
Resulting map looks similar, only shading has moved a bit. I'm not sure
if this is the result of differences in interpolation or maybe change
introduced by moving DEM borders.
--
Best regards,
Andrzej
Index: src/uk/me/parabola/imgfmt/app/dem/DEMFile.java
===
--- src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (revision 4056)
+++ src/uk/me/parabola/imgfmt/app/dem/DEMFile.java (working copy)
@@ -22,6 +22,7 @@
import uk.me.parabola.imgfmt.app.ImgFileWriter;
import uk.me.parabola.imgfmt.fs.ImgChannel;
import uk.me.parabola.mkgmap.reader.hgt.HGTConverter;
+import uk.me.parabola.imgfmt.Utils;
/**
* The DEM file. This consists of information about elevation. It is used for
hill shading
@@ -55,6 +56,26 @@
HGTConverter hgtConverter = new HGTConverter(pathToHGT, area,
demPolygonMapUnits);
hgtConverter.setOutsidePolygonHeight(outsidePolygonHeight);
+ int res = hgtConverter.getHighestRes();
+ if (res > 0) {
+ if (pointDistances.get(0) == -1 ||
pointDistances.get(0) == (int) ((1 << 29) / (res * 45))) {
+ // align area with HGT raster
+ double hgtDis = 1.0D/res;
+ double top = Utils.toDegrees(area.getMaxLat());
+ double bottom =
Utils.toDegrees(area.getMinLat());
+ double left =
Utils.toDegrees(area.getMinLong());
+ double right =
Utils.toDegrees(area.getMaxLong());
+
+ top = Math.ceil(top/hgtDis)*hgtDis;
+ bottom = Math.floor(bottom/hgtDis)*hgtDis;
+ left = Math.floor(left/hgtDis)*hgtDis;
+ right = Math.ceil(right/hgtDis)*hgtDis;
+
+ area = new Area(bottom, left, top, right);
+ //System.out.println("DEM aligned to HGT");
+ }
+ }
+
int zoom = 0;
int lastDist = pointDistances.get(pointDistances.size()-1);
for (int pointDist : pointDistances) {
___
mkgmap-dev mailing list
mkgmap-dev@lists.mkgmap.org.uk
http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, yes, that looks strange. I've added option --show-profiles=1 to yours so that MapSource displays height values even in those areas where hill shading is not working. For all: The problem is near N49.19004 E20.07117 Mapsource displays no hill shading at zoom 500m, the wrong area changes when you zoom in and out or press 2xStrg+G to clear the cache. Not sure if this is a problem in the DEM data or in the hill shading algo. Will try to find out more tomorrow. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Donnerstag, 11. Januar 2018 16:50:13 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I have uploaded at http://files.mkgmap.org.uk/ data for compilation, which gives me a map with shading errors. Mkgmap version 4045. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I have uploaded at http://files.mkgmap.org.uk/ data for compilation, which gives me a map with shading errors. Mkgmap version 4045. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I will try to prepare a simple example. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, please let us know how to reproduce the problem maybe post the problem tile. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Dienstag, 9. Januar 2018 22:11:12 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, yes, the problem is still present in 4041. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Am 10.01.2018 um 18:15 schrieb Gerd Petermann: Hi Peter, I think the only alternative to interpolation is to use the height of the nearest hgt point for a given DEM point. That is my theory, but I can't prove lg Peter One may minimize this distance by chosing proper dem-dist values and by aligning the tiles to the raster given by hgt. Example: If the upper left corner of your img file is at lat 51.0 , lon 10.0 and dem-dist is 2^32/(360*3600) ~ 3314, the DEM points will be close to the hgt points in file N51E010.hgt if that is in 1'' resolution. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Peter Danninger <pe...@danninger.eu> Gesendet: Mittwoch, 10. Januar 2018 17:53:46 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() I don't think rounding or interpolating will deliver any better results. If you have 1" .hgt files, you get an altitude value every 30m. With 3" .hgt files much worse, every 90m. In the mountains, between 2 points may be a big wall with some 100m height-difference. So if you have small deviations in lat/lon, you will evtl. get very false elevation values. If you interpolate altitude-values, you calculate wrong values by design :-( May be I'm wrong, but this is my theory. lg Peter Am 10.01.2018 um 17:39 schrieb Gerd Petermann: Hi Andrzej, if you want to try this: you can easily change the code in interpolatedHeight to return the height in feet, just make sure that you don't convert UNDEF. The only other change that is needed is in this line in DemHeader: writer.putInt(0); // 0: elevation in metres, 1: foot well, the comment should say feet , not foot Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Frank Stinner <frank.stin...@kabelmail.de> Gesendet: Mittwoch, 10. Januar 2018 17:27:22 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd and Andrzej, i think "overkill" is a good word. For algorithm that'a only numbers, it does not matter. The questions is, how exact are the hgt-values. We don't know that, but i don't believe it is +-1m or so. I'm not wondering, when the hgt's have +-5m or +-10m. The copernicus-data have +-7m and i don't believe the technic was worse. That's why feets are overkill. The numbers are 3 times greater, that's why the dem's are greater. It's not worth it. Frank --- Diese E-Mail wurde von Avast Antivirus-Software auf Viren geprüft. https://www.avast.com/antivirus ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi, > If you interpolate altitude-values, you calculate wrong values > by design :-( Would it be possible, to create first layer without interpolation? On my maps I'm trying to set DEM resolution equal to HGT, but I can't control offset and I'm not sure if settings are precise enough to maintain pixel to pixel accuracy. I would like to get a mode, where HGT are encoded without resampling. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Peter, I think the only alternative to interpolation is to use the height of the nearest hgt point for a given DEM point. One may minimize this distance by chosing proper dem-dist values and by aligning the tiles to the raster given by hgt. Example: If the upper left corner of your img file is at lat 51.0 , lon 10.0 and dem-dist is 2^32/(360*3600) ~ 3314, the DEM points will be close to the hgt points in file N51E010.hgt if that is in 1'' resolution. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Peter Danninger <pe...@danninger.eu> Gesendet: Mittwoch, 10. Januar 2018 17:53:46 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() I don't think rounding or interpolating will deliver any better results. If you have 1" .hgt files, you get an altitude value every 30m. With 3" .hgt files much worse, every 90m. In the mountains, between 2 points may be a big wall with some 100m height-difference. So if you have small deviations in lat/lon, you will evtl. get very false elevation values. If you interpolate altitude-values, you calculate wrong values by design :-( May be I'm wrong, but this is my theory. lg Peter Am 10.01.2018 um 17:39 schrieb Gerd Petermann: > Hi Andrzej, > > if you want to try this: > you can easily change the code in interpolatedHeight to return the height in > feet, just make sure that you don't convert UNDEF. > > The only other change that is needed is in this line in DemHeader: > writer.putInt(0); // 0: elevation in metres, 1: foot > > well, the comment should say feet , not foot > > Gerd > > > Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Frank > Stinner <frank.stin...@kabelmail.de> > Gesendet: Mittwoch, 10. Januar 2018 17:27:22 > An: mkgmap-dev@lists.mkgmap.org.uk > Betreff: Re: [mkgmap-dev] HGT - getElevation() > > Hi Gerd and Andrzej, > > i think "overkill" is a good word. > > For algorithm that'a only numbers, it does not matter. The questions is, how > exact are the hgt-values. We don't know that, but i don't believe it is > +-1m or so. I'm not wondering, when the hgt's have +-5m or +-10m. The > copernicus-data have +-7m and i don't believe the technic was worse. > > That's why feets are overkill. The numbers are 3 times greater, that's why > the dem's are greater. It's not worth it. > > > Frank > > --- > Diese E-Mail wurde von Avast Antivirus-Software auf Viren geprüft. > https://www.avast.com/antivirus > > > ___ > mkgmap-dev mailing list > mkgmap-dev@lists.mkgmap.org.uk > http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev > ___ > mkgmap-dev mailing list > mkgmap-dev@lists.mkgmap.org.uk > http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev > ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
I don't think rounding or interpolating will deliver any better results. If you have 1" .hgt files, you get an altitude value every 30m. With 3" .hgt files much worse, every 90m. In the mountains, between 2 points may be a big wall with some 100m height-difference. So if you have small deviations in lat/lon, you will evtl. get very false elevation values. If you interpolate altitude-values, you calculate wrong values by design :-( May be I'm wrong, but this is my theory. lg Peter Am 10.01.2018 um 17:39 schrieb Gerd Petermann: Hi Andrzej, if you want to try this: you can easily change the code in interpolatedHeight to return the height in feet, just make sure that you don't convert UNDEF. The only other change that is needed is in this line in DemHeader: writer.putInt(0); // 0: elevation in metres, 1: foot well, the comment should say feet , not foot Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Frank Stinner <frank.stin...@kabelmail.de> Gesendet: Mittwoch, 10. Januar 2018 17:27:22 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd and Andrzej, i think "overkill" is a good word. For algorithm that'a only numbers, it does not matter. The questions is, how exact are the hgt-values. We don't know that, but i don't believe it is +-1m or so. I'm not wondering, when the hgt's have +-5m or +-10m. The copernicus-data have +-7m and i don't believe the technic was worse. That's why feets are overkill. The numbers are 3 times greater, that's why the dem's are greater. It's not worth it. Frank --- Diese E-Mail wurde von Avast Antivirus-Software auf Viren geprüft. https://www.avast.com/antivirus ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, if you want to try this: you can easily change the code in interpolatedHeight to return the height in feet, just make sure that you don't convert UNDEF. The only other change that is needed is in this line in DemHeader: writer.putInt(0); // 0: elevation in metres, 1: foot well, the comment should say feet , not foot Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Frank Stinner <frank.stin...@kabelmail.de> Gesendet: Mittwoch, 10. Januar 2018 17:27:22 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd and Andrzej, i think "overkill" is a good word. For algorithm that'a only numbers, it does not matter. The questions is, how exact are the hgt-values. We don't know that, but i don't believe it is +-1m or so. I'm not wondering, when the hgt's have +-5m or +-10m. The copernicus-data have +-7m and i don't believe the technic was worse. That's why feets are overkill. The numbers are 3 times greater, that's why the dem's are greater. It's not worth it. Frank --- Diese E-Mail wurde von Avast Antivirus-Software auf Viren geprüft. https://www.avast.com/antivirus ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd and Andrzej, i think "overkill" is a good word. For algorithm that'a only numbers, it does not matter. The questions is, how exact are the hgt-values. We don't know that, but i don't believe it is +-1m or so. I'm not wondering, when the hgt's have +-5m or +-10m. The copernicus-data have +-7m and i don't believe the technic was worse. That's why feets are overkill. The numbers are 3 times greater, that's why the dem's are greater. It's not worth it. Frank --- Diese E-Mail wurde von Avast Antivirus-Software auf Viren geprüft. https://www.avast.com/antivirus ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, height precision is 3 times better in feet and rounding errors smaller. It is probably overkill, but might be handy with good data. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, I can easily change the code to support feet. I think Franks program already supports that. But what would be the benefit ? Even with 1'' hgt data we have only one data point every ~25 metres (in western europe) giving the height in metres, I can't believe that we would improve quality by using feet, but runtime and DEM file size are both much higher. I've justed tested it with my DEMCompressor: Compress 32 1'' files at N60 with metre takes 28 sec and creates ~65M *.DEM files. Same files with feet requires 32 secs and creates ~105 M *.dem files. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Mittwoch, 10. Januar 2018 16:14:15 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, no idea about int vs. short, but I wouldn't expect any differences on Intel CPU. I think it could be interesting to return interpolated values in feet, assuming that you can encode whole DEM in feet. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, no idea about int vs. short, but I wouldn't expect any differences on Intel CPU. I think it could be interesting to return interpolated values in feet, assuming that you can encode whole DEM in feet. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, sorry for my stupid error and thanks for the patch. I kept the double because I compared the results with those from class InterpolationBilinear, but I agree this should not be used in the final code. See http://www.mkgmap.org.uk/websvn/revision.php?repname=mkgmap=4043 Maybe I should change all values to int, maybe that is a bit faster. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Mittwoch, 10. Januar 2018 14:33:40 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd, I see you have introduced NO_VAL constant in recent code. Please note, then my interpolation code can return (double)HGTReader.UNDEF too, when it returns nearest HGT value, which happened to be void. I think more clear would be to move Math.round() to interpolatedHeight(), see attached patch. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
forgotten patch here.
--
Best regards,
Andrzej
Index: src/uk/me/parabola/mkgmap/reader/hgt/HGTConverter.java
===
--- src/uk/me/parabola/mkgmap/reader/hgt/HGTConverter.java (revision 4042)
+++ src/uk/me/parabola/mkgmap/reader/hgt/HGTConverter.java (working copy)
@@ -112,8 +112,7 @@
int hRB = rdr.ele(xRight, yBottom);
lastRow = row;
- double rc = interpolatedHeight(x1-xLeft, y1-yBottom, hLT, hRT,
hRB, hLB);
- short rc2 = (rc != NO_VAL) ? ((short) Math.round(rc)) :
HGTReader.UNDEF;
+ short rc2 = interpolatedHeight(x1-xLeft, y1-yBottom, hLT, hRT,
hRB, hLB);
// double lon = lon32 * FACTOR;
// double lat = lat32 * FACTOR;
// System.out.println(String.format("%.7f %.7f: (%.2f) %d",
lat,lon,rc,rc2));
@@ -149,60 +148,60 @@
* @param hlb height left bottom
* @return the interpolated height
*/
- private static double interpolatedHeight(double qx, double qy, int hlt,
int hrt, int hrb, int hlb) {
+ private static short interpolatedHeight(double qx, double qy, int hlt,
int hrt, int hrb, int hlb) {
// extrapolate single node height if requested point is not near
// for multiple missing nodes, return the height of the neares
node
if (hlb == HGTReader.UNDEF) {
if (hrb == HGTReader.UNDEF || hlt == HGTReader.UNDEF ||
hrt == HGTReader.UNDEF) {
if (hrt != HGTReader.UNDEF && hlt !=
HGTReader.UNDEF && qy > 0.5D) //top edge
- return (1.0D - qx)*hlt + qx*hrt;
+ return (short)Math.round((1.0D -
qx)*hlt + qx*hrt);
if (hrt != HGTReader.UNDEF && hrb !=
HGTReader.UNDEF && qx > 0.5D) //right edge
- return (1.0D - qy)*hrb + qy*hrt;
+ return (short)Math.round((1.0D -
qy)*hrb + qy*hrt);
//if (hlt != HGTReader.UNDEF && hrb !=
HGTReader.UNDEF && qx + qy > 0.5D && gx + qy < 1.5D) //diagonal
// nearest value
- return (double)((qx < 0.5D)? ((qy < 0.5D)? hlb:
hlt): ((qy < 0.5D)? hrb: hrt));
+ return (short)((qx < 0.5D)? ((qy < 0.5D)? hlb:
hlt): ((qy < 0.5D)? hrb: hrt));
}
if (qx + qy < 0.4D) // point is near missing value
- return NO_VAL;
+ return HGTReader.UNDEF;
hlb = hlt + hrb - hrt;
} else if (hrt == HGTReader.UNDEF) {
if (hlb == HGTReader.UNDEF || hrb == HGTReader.UNDEF ||
hlt == HGTReader.UNDEF) {
if (hlb != HGTReader.UNDEF && hrb !=
HGTReader.UNDEF && qy < 0.5D) //lower edge
- return (1.0D - qx)*hlb + qx*hrb;
+ return (short)Math.round((1.0D -
qx)*hlb + qx*hrb);
if (hlb != HGTReader.UNDEF && hlt !=
HGTReader.UNDEF && qx < 0.5D) //left edge
- return (1.0D - qy)*hlb + qy*hlt;
+ return (short)Math.round((1.0D -
qy)*hlb + qy*hlt);
//if (hlt != HGTReader.UNDEF && hrb !=
HGTReader.UNDEF && qx + qy > 0.5D && gx + qy < 1.5D) //diagonal
// nearest value
- return (double)((qx < 0.5D)? ((qy < 0.5D)? hlb:
hlt): ((qy < 0.5D)? hrb: hrt));
+ return (short)((qx < 0.5D)? ((qy < 0.5D)? hlb:
hlt): ((qy < 0.5D)? hrb: hrt));
}
if (qx + qy > 1.6D) // point is near missing value
- return NO_VAL;
+ return HGTReader.UNDEF;
hrt = hlt + hrb - hlb;
} else if (hrb == HGTReader.UNDEF) {
if (hlb == HGTReader.UNDEF || hlt == HGTReader.UNDEF ||
hrt == HGTReader.UNDEF) {
if (hlt != HGTReader.UNDEF && hrt !=
HGTReader.UNDEF && qy > 0.5D) //top edge
- return (1.0D - qx)*hlt + qx*hrt;
+ return (short)Math.round((1.0D -
qx)*hlt + qx*hrt);
if (hlt != HGTReader.UNDEF && hlb !=
HGTReader.UNDEF && qx < 0.5D) //left edge
- return (1.0D - qy)*hlb + qy*hlt;
+ return (short)Math.round((1.0D -
qy)*hlb + qy*hlt);
//if (hlb != HGTReader.UNDEF &&
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, I see you have introduced NO_VAL constant in recent code. Please note, then my interpolation code can return (double)HGTReader.UNDEF too, when it returns nearest HGT value, which happened to be void. I think more clear would be to move Math.round() to interpolatedHeight(), see attached patch. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Henning, hgt files contain voids, and we don't write 0 height for such a void, instead we write a value that is considered to be invalid. When you hover over such a point MapSource doesn't display a height. I did not yet try to find out how this influences the ele values in calculated gpx tracks. I assume that a gpx point without ele attribute is created. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Henning Scholland <o...@hscholland.de> Gesendet: Mittwoch, 10. Januar 2018 13:48:11 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() On 10.01.2018 20:33, Andrzej Popowski wrote: > > i would be very defensive with interpolation in the case of > > missing values. > > I got the same feeling. The best way to fill voids is to process whole > HGT. And I believe that DEM providers already have done it, so any > attempt at simple extrapolation would rather add errors than make > output better. In my code I preserve voids. If requested coordinates > are near void HGT node, then returned value is void. > Hi Andrzej, >From user perspective, I think it's more correct to somehow fill voids based of surrounding data than using elevation = 0. Elevation=0 I would only use as a fall back, if whole hgt-file is missing. I agree with you, that we shouldn't spent a lot of code for guessing. To have a simple example: If there would be a cone-shape mountain in reality, but in DEM the top of the cone is missing, I would suggest to fill the hole with same elevation as the highest circle. So result in DEM would be a truncated cone shape. Of course reality is not that simple. :'( ;-) Henning ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
On 10.01.2018 20:33, Andrzej Popowski wrote: > > i would be very defensive with interpolation in the case of > > missing values. > > I got the same feeling. The best way to fill voids is to process whole > HGT. And I believe that DEM providers already have done it, so any > attempt at simple extrapolation would rather add errors than make > output better. In my code I preserve voids. If requested coordinates > are near void HGT node, then returned value is void. > Hi Andrzej, >From user perspective, I think it's more correct to somehow fill voids based of surrounding data than using elevation = 0. Elevation=0 I would only use as a fall back, if whole hgt-file is missing. I agree with you, that we shouldn't spent a lot of code for guessing. To have a simple example: If there would be a cone-shape mountain in reality, but in DEM the top of the cone is missing, I would suggest to fill the hole with same elevation as the highest circle. So result in DEM would be a truncated cone shape. Of course reality is not that simple. :'( ;-) Henning ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Frank and Gerd, > It's a little bit more to calculate The actual code is quite simple, it is 3 times linear interpolation: double hxt = (1.0D - qx)*hlt + qx*hrt; double hxb = (1.0D - qx)*hlb + qx*hrb; return (1.0D - qy)*hxb + qy*hxt; > i would be very defensive with interpolation in the case of > missing values. I got the same feeling. The best way to fill voids is to process whole HGT. And I believe that DEM providers already have done it, so any attempt at simple extrapolation would rather add errors than make output better. In my code I preserve voids. If requested coordinates are near void HGT node, then returned value is void. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Frank what problem do you see when we interpolate data? What is the advantage of voids? Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Frank Stinner <frank.stin...@leipzig.de> Gesendet: Mittwoch, 10. Januar 2018 10:53:44 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi Gerd and Andrzej, > your choice of triangle is arbitrary. of course, but the version with 4 triangles isn't. > bilinear interpolation seems better Yes. It's a little bit more to calculate, but better. For our special case it is: h = h11 * (1 + y * x - x - y) + h21 * (x - y * x) + h12 * (y - y * x) + h22 * y * x (h11 is left-bottom in coordinate origin, the horizontal and vertical distance is 1) If your java-class have to much overhead, you can use this formula. By the way, i would be very defensive with interpolation in the case of missing values. The hgt values are only interpolated values from original measurments. I would not interpolate values in the near of missing values. That give us a nice but unreliably picture. Frank ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd and Andrzej, > your choice of triangle is arbitrary. of course, but the version with 4 triangles isn't. > bilinear interpolation seems better Yes. It's a little bit more to calculate, but better. For our special case it is: h = h11 * (1 + y * x - x - y) + h21 * (x - y * x) + h12 * (y - y * x) + h22 * y * x (h11 is left-bottom in coordinate origin, the horizontal and vertical distance is 1) If your java-class have to much overhead, you can use this formula. By the way, i would be very defensive with interpolation in the case of missing values. The hgt values are only interpolated values from original measurments. I would not interpolate values in the near of missing values. That give us a nice but unreliably picture. Frank___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd, yes, the problem is still present in 4041. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, the patch is based on r4040. Do you see those artifacts with r4041? Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Dienstag, 9. Januar 2018 21:04:20 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi, Gerd, your additional interpolation would be redundant. Side note: comparison of double and integer doesn't look safe. Maybe would be better if interpolateHeight() returned integer? Frank, I understand your idea, but please note, that your choice of triangle is arbitrary. There are 2 ways of dividing a square and if you would choose the other way, results would be different and not only for the middle point. This suggest that the solution isn't complete. Classic bilinear interpolation seems better, but differences are small and actual map looks more or less the same. I have attached second patch, which apply interpolation when only 2 values on any edge are present. Another topic. I observe some DEM artifacts, when displaying area near the border of dem-poly. There are small rectangles without shading within a bigger tile. They appear at random, when scrolling or zooming map in BaseCamp or Mapsource. Ctrl-G doesn't help. See attached img. I haven't noticed the problem, when compiling a map without dem-poly. Maybe it is a result of slanted edge of DEM clipped by poly. -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi,
Gerd, your additional interpolation would be redundant.
Side note: comparison of double and integer doesn't look safe. Maybe
would be better if interpolateHeight() returned integer?
Frank, I understand your idea, but please note, that your choice of
triangle is arbitrary. There are 2 ways of dividing a square and if you
would choose the other way, results would be different and not only for
the middle point. This suggest that the solution isn't complete. Classic
bilinear interpolation seems better, but differences are small and
actual map looks more or less the same.
I have attached second patch, which apply interpolation when only 2
values on any edge are present.
Another topic. I observe some DEM artifacts, when displaying area near
the border of dem-poly. There are small rectangles without shading
within a bigger tile. They appear at random, when scrolling or zooming
map in BaseCamp or Mapsource. Ctrl-G doesn't help. See attached img.
I haven't noticed the problem, when compiling a map without dem-poly.
Maybe it is a result of slanted edge of DEM clipped by poly.
--
Best regards,
Andrzej
Index: src/uk/me/parabola/mkgmap/reader/hgt/HGTConverter.java
===
--- src/uk/me/parabola/mkgmap/reader/hgt/HGTConverter.java (revision 4040)
+++ src/uk/me/parabola/mkgmap/reader/hgt/HGTConverter.java (working copy)
@@ -109,7 +109,7 @@
int hLB = rdr.ele(xLeft, yBottom);
int hRB = rdr.ele(xRight, yBottom);
lastRow = row;
- double rc = interpolatedHeightInNormatedRectangle(x1-xLeft,
y1-yBottom, hLT, hRT, hRB, hLB);
+ double rc = interpolatedHeight(x1-xLeft, y1-yBottom, hLT, hRT,
hRB, hLB);
if (rc == HGTReader.UNDEF) {
int sum = 0;
int valid = 0;
@@ -147,10 +147,10 @@
}
}
}
-
+
/**
* Interpolate the height of point p from the 4 closest values in the
hgt matrix.
-* Code is a copy from Frank Stinners program BuildDEMFile
(Hgtreader.cs)
+* Bilinear interpolation with single node restore
* @param qx value from 0 .. 1 gives relative x position in matrix
* @param qy value from 0 .. 1 gives relative y position in matrix
* @param hlt height left top
@@ -159,42 +159,70 @@
* @param hlb height left bottom
* @return the interpolated height
*/
- double interpolatedHeightInNormatedRectangle(double qx, double qy, int
hlt, int hrt, int hrb, int hlb) {
- if (hlb == HGTReader.UNDEF || hrt == HGTReader.UNDEF)
- return HGTReader.UNDEF; // keine Berechnung mУЖglich
-
- /*
-* In welchem Dreieck liegt der Punkt? oben +-/ |/
-*
-* unten /| /-+
-*/
- if (qy >= qx) { // oberes Dreieck aus hlb, hrt und hlt (Anstieg
py/px
- // ist grУЖУer als
height/width)
-
- if (hlt == HGTReader.UNDEF)
+ double interpolatedHeight(double qx, double qy, int hlt, int hrt, int
hrb, int hlb) {
+ // extrapolate single node height if requested point is not near
+ // for multiple missing nodes, return the height of the neares
node
+ if (hlb == HGTReader.UNDEF) {
+ if (hrb == HGTReader.UNDEF || hlt == HGTReader.UNDEF ||
hrt == HGTReader.UNDEF) {
+ if (hrt != HGTReader.UNDEF && hlt !=
HGTReader.UNDEF && qy > 0.5D) //top edge
+ return (1.0D - qx)*hlt + qx*hrt;
+ if (hrt != HGTReader.UNDEF && hrb !=
HGTReader.UNDEF && qx > 0.5D) //right edge
+ return (1.0D - qy)*hrb + qy*hrt;
+ //if (hlt != HGTReader.UNDEF && hrb !=
HGTReader.UNDEF && qx + qy > 0.5D && gx + qy < 1.5D) //diagonal
+ // nearest value
+ return (double)((qx < 0.5D)? ((qy < 0.5D)? hlb:
hlt): ((qy < 0.5D)? hrb: hrt));
+ }
+ if (qx + qy < 0.4D) // point is near missing value
return HGTReader.UNDEF;
-
- // hlt als Koordinatenursprung normieren; mit hrt und
hlb 3 Punkte
- // einer Ebene (3-Punkt-Gleichung)
- hrt -= hlt;
- hlb -= hlt;
- qy -= 1;
-
- return hlt + qx * hrt - qy * hlb;
-
- } else { // unteres Dreieck aus hlb, hrb und hrt
-
- if (hrb == HGTReader.UNDEF)
+ hlb = hlt + hrb - hrt;
+ }
+ else if (hrt == HGTReader.UNDEF) {
+
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej and Gerd, for your examples 1 0 0 7 0 1 7 0 is 2 for the middlepoint mathematical correct: (h1 + h2 + h3 + h4) / 4. My interpolation is the that: We live in a world with "triangular" surface (like in a computergame). We looking for a point p for the "surrounding" 4 points from the hgt. They form a rectangle (or square). The 4 point form 2 triangles (for me): hlthrt (height right top) +---+ | /| | / | |/ | +---+ hlbhrb We looking in which triangle our point p is. A triangle define a plane and we can use a "3-Punkt-Gleichung" (don't know the english word). For the left triangle: use hlt as coordinate origin hrt -= hlt; hlb -= hlt; qy -= 1; and calculate hlt + qx * hrt - qy * hlb For the right triangle: use hrb as coordinate origin hrt -= hrb; hlb -= hrb; qx -= 1; and calculate hrb - qx * hlb + qy * hrt This principle can be a little extend: hlthrt +---+ |\ /| | x | |/ \| +---+ hlbhrb It's easy to calculate the middlepoint and then we have 4 triangles. Then we have to decide, which triangle we need. Just we have 4 triangles. (So we have our "triangular" surface.) And so on ... Frank --- Diese E-Mail wurde von Avast Antivirus-Software auf Viren geprüft. https://www.avast.com/antivirus ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej,
the patch keeps my simple code for an arithmetic mean if interpolatedHeight
returns UNDEF:
if (rc == HGTReader.UNDEF) {
int sum = 0;
int valid = 0;
int[] heights = { hLT, hRT, hLB, hRB };
for (int h : heights) {
if (h == HGTReader.UNDEF)
continue;
valid++;
sum += h;
}
if(valid >= 2)
rc = Math.round((double)sum/valid);
}
I think this should be removed?
Gerd
Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej
Popowski <po...@poczta.onet.pl>
Gesendet: Dienstag, 9. Januar 2018 19:10:15
An: mkgmap-dev@lists.mkgmap.org.uk
Betreff: Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
I have prepared a patch with bilinear interpolation, which supports
missing values. It can restore single missing value, but if more values
are missing it simply returns the nearest value without interpolation.
--
Best regards,
Andrzej
___
mkgmap-dev mailing list
mkgmap-dev@lists.mkgmap.org.uk
http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Andrzej, Frank,
I found class InterpolationBilinear
http://download.java.net/media/jai/javadoc/1.1.3/jai-apidocs/javax/media/jai/InterpolationBilinear.html
which seems to do exactly what you propose when all four values are valid.
I found no code for triangular interpolation yet, but I assume we can use the
existing code in interpolatedHeightInNormatedRectangle with some more math.
We just have to "rotate" the values and xchange qx and qy so that
interpolatedHeightInNormatedRectangle never returns UNDEF if three values are
valid.
Right?
Gerd
Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej
Popowski <po...@poczta.onet.pl>
Gesendet: Dienstag, 9. Januar 2018 15:45:14
An: mkgmap-dev@lists.mkgmap.org.uk
Betreff: Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
I'm trying to find if resolution can be better. I think current
interpolation is not good. Imagine two combination of 4-points heights
from HGT with values like:
1 0
0 7
and
0 1
7 0
For a point in the middle of rectangle, current interpolation would give
result 0 in first case and 4 in second, while I would expect 2 in both
cases.
I think the code, which calculates interpolation based on 3 points,
should be used only if the forth is void. There could be 3 combinations
of triangles in this case. But if all 4 points are valid, there should
be applied standard bilinear interpolation. Something like:
if (hlb == HGTReader.UNDEF) {
return triangle(...
}
if (hrb == HGTReader.UNDEF) {
return triangle(...
}
if (hlt == HGTReader.UNDEF) {
return triangle(...
}
if (hrt == HGTReader.UNDEF) {
return triangle(...
}
// 4 points valid
return bilinear(...
--
Best regards,
Andrzej
___
mkgmap-dev mailing list
mkgmap-dev@lists.mkgmap.org.uk
http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
___
mkgmap-dev mailing list
mkgmap-dev@lists.mkgmap.org.uk
http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej,
to be honest I did not understand the math, that's why I've copied the code
from Frank.
Maybe I've done something wrong with the calculation of dx and dy.
Do you think that Franks code has the same problem?
Gerd
Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej
Popowski <po...@poczta.onet.pl>
Gesendet: Dienstag, 9. Januar 2018 15:45:14
An: mkgmap-dev@lists.mkgmap.org.uk
Betreff: Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
I'm trying to find if resolution can be better. I think current
interpolation is not good. Imagine two combination of 4-points heights
from HGT with values like:
1 0
0 7
and
0 1
7 0
For a point in the middle of rectangle, current interpolation would give
result 0 in first case and 4 in second, while I would expect 2 in both
cases.
I think the code, which calculates interpolation based on 3 points,
should be used only if the forth is void. There could be 3 combinations
of triangles in this case. But if all 4 points are valid, there should
be applied standard bilinear interpolation. Something like:
if (hlb == HGTReader.UNDEF) {
return triangle(...
}
if (hrb == HGTReader.UNDEF) {
return triangle(...
}
if (hlt == HGTReader.UNDEF) {
return triangle(...
}
if (hrt == HGTReader.UNDEF) {
return triangle(...
}
// 4 points valid
return bilinear(...
--
Best regards,
Andrzej
___
mkgmap-dev mailing list
mkgmap-dev@lists.mkgmap.org.uk
http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
___
mkgmap-dev mailing list
mkgmap-dev@lists.mkgmap.org.uk
http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Gerd,
I'm trying to find if resolution can be better. I think current
interpolation is not good. Imagine two combination of 4-points heights
from HGT with values like:
1 0
0 7
and
0 1
7 0
For a point in the middle of rectangle, current interpolation would give
result 0 in first case and 4 in second, while I would expect 2 in both
cases.
I think the code, which calculates interpolation based on 3 points,
should be used only if the forth is void. There could be 3 combinations
of triangles in this case. But if all 4 points are valid, there should
be applied standard bilinear interpolation. Something like:
if (hlb == HGTReader.UNDEF) {
return triangle(...
}
if (hrb == HGTReader.UNDEF) {
return triangle(...
}
if (hlt == HGTReader.UNDEF) {
return triangle(...
}
if (hrt == HGTReader.UNDEF) {
return triangle(...
}
// 4 points valid
return bilinear(...
--
Best regards,
Andrzej
___
mkgmap-dev mailing list
mkgmap-dev@lists.mkgmap.org.uk
http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi Andrzej, thanks for reviewing the code. Yes, they should never be negative. Gerd Von: mkgmap-dev <mkgmap-dev-boun...@lists.mkgmap.org.uk> im Auftrag von Andrzej Popowski <po...@poczta.onet.pl> Gesendet: Dienstag, 9. Januar 2018 13:30:35 An: mkgmap-dev@lists.mkgmap.org.uk Betreff: Re: [mkgmap-dev] HGT - getElevation() Hi, correction - these values never are negative? -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
Re: [mkgmap-dev] HGT - getElevation()
Hi, correction - these values never are negative? -- Best regards, Andrzej ___ mkgmap-dev mailing list mkgmap-dev@lists.mkgmap.org.uk http://www.mkgmap.org.uk/mailman/listinfo/mkgmap-dev
