Gigas002 / GTiff2Tiles

Analogue of gdal2tiles.py/MapTiler on C#.
https://gigas002.github.io/GTiff2Tiles/
Do What The F*ck You Want To Public License
15 stars 4 forks source link

TileGrid Row issue after zoom 8 #9

Closed carnegiea closed 5 years ago

carnegiea commented 5 years ago

The tilegrid after zoom 8 for my image (tiff) which is located near to equator is not correct when refering to the grid at https://tile.gbif.org/ui/4326/.

Lets say for the X: 404, Y:127 at Z:8 The website reference showing 127 but from the TilesMinMax in SetCropProperties method in Image.cs file shows X:404, Y:129 at Z: 8. And it gets increasingly large after that.

I have added these lines in Tile.cs

if (zoom < 8) { tilesYs[0] = Convert.ToInt32(Math.Ceiling((90.0 + minY) Math.Pow(2.0, zoom) / 180.0) - 1.0); tilesYs[1] = Convert.ToInt32(Math.Ceiling((90.0 + maxY) Math.Pow(2.0, zoom) / 180.0) - 1.0); } else { tilesYs[0] = Convert.ToInt32(Math.Ceiling((90.0 + minY) Math.Pow(2.0, zoom) / 180.0) - 1.0 - (Math.Pow(2.0, 2 + (zoom - 8)) - 2)); tilesYs[1] = Convert.ToInt32(Math.Ceiling((90.0 + maxY) Math.Pow(2.0, zoom) / 180.0) - 1.0 - (Math.Pow(2.0, 2 + (zoom - 8)) - 2));

            }

and found out the output same as the one from the https://tile.gbif.org/ui/4326/.

But... later on errors are triggered at GeoQuery in Image.cs where after zoom 8 the readPosMinY will reach the maximum Y and therefore lead to readYSize become 0, and trigger error at

NetVips.Image tileImage; try { tileImage = inputImage.Crop(readPosX, readPosY, readXSize, readYSize); } catch (Exception exception) { throw new ImageException(string.Format(Strings.UnableToCreateTile, tileX, tileY), exception); }

which is at line 238 in Image.cs

Kindly advise if I miss anything or misunderstood anything regarding this matter.

Thanks a lot.

Gigas002 commented 5 years ago

Hello. That issue is probably happens because that website (I'm not sure) doesn't use tms tile naming system (even though they are EPSG:4326). The difference is indeed in y tile numbers. GTiff2Tiles now supports only tms naming system. Here you can read more about it. To explicitly flip y number from code, you can use this:

//Tile.GetTileNumbersFromCoords method, 
int y1 = Convert.ToInt32(Math.Pow(2.0, zoom) - tilesYs.Min() - 1);
int y2 = Convert.ToInt32(Math.Pow(2.0, zoom) - tilesYs.Max() - 1);
int[] ys = {y1, y2};
return (tilesXs.Min(), ys.Min(), tilesXs.Max(), ys.Max());

Then, you should also change isFlipY parameter in Image.WriteTile (the one with NetVips.Image parameter) method:

//226 string in Image.cs
(double minX, double minY, double maxX, double maxY) = Tile.Tile.TileBounds(tileX, tileY, zoom, true);

It should work for crop algorithm (didn't test on join yet).

Hope this helps you. Please, let me know if it solves your issue, and also feel free to ask more questions if you need. I think I will add functionality for non-tms tiles in the future.

carnegiea commented 5 years ago

The web reference https://tile.gbif.org/ui/4326/ is similar to what is implemented in Openlayers and I am not sure if openlayer also implementing non-tms naming system.

By the way the above modifications work but now the Y naming is falling into different X folders.

Using the results as similar to modifications above, the layer could not be overlaid as the same location as the base in openlayers.

Please let me know if you need a tiff file to test it out?

Also, the output png is not showing any content, unless i change to jpegsave.

Thanks

Gigas002 commented 5 years ago

Yes, that would be great, if you're able to send me that .tif. I tested the above code on my data, and the output is the same, as MapTiler Pro creates (with the following params: -geodetic -o D:/Test/MaptilerResult -work_dir D:/Test/MaptilerTemp -srs "+proj=longlat +datum=WGS84 +no_defs" -zoom 0 11 D:/Test/map_16_4326.tif (throwing out the -tms argument)).

P.S. I found this code in gdal2tiles.py:

if self.options.tmscompatible:
    args['tmsoffset'] = "-1"
else:
    args['tmsoffset'] = ""

I suppose, there's exactly 1 tile number for y difference between tms tiles output only for OpenLayers. So, changing outputTileFileInfo value in Image.WriteTile method like this:

FileInfo outputTileFileInfo = new FileInfo(Path.Combine(tileDirectoryInfo.FullName, $"{tileY + 1}{Enums.Extensions.Png}"));

Gives me the same result, as on that website.

carnegiea commented 5 years ago

The link to the tiff is at (~720MB)

https://1drv.ms/u/s!AmpFsfFEKx7Ki4Va0rdy39mhQYCXCw

Perhaps you can test the tiff above with openlayer if possible. For this particular tiff location, the tileY is deviating from the openlayer grid by a series of increment value.

I will download a new set of GTiff2Tiles and do test again on this matter, then will have a clearer picture on what is happening.

Thanks a lot for your help!

carnegiea commented 5 years ago

Refering to this https://www.maptiler.com/google-maps-coordinates-tile-bounds-projection/

Zoom 8 , X: 202, Y:127

I have downloaded a new set and run both test on TMS turned on and off, the X folder seems incorrect.

RUN ON CROPPING:

I put the maptiler result as reference here Google Zoom 8: X: 202, Y: 127 Zoom 9: X: 404, Y: 254 Zoom 10: X: 808, Y: 508 Zoom 11: X: 1617, Y: 1016 Zoom 12: X: 3234, Y: 2031

TMS Zoom 8: X: 202, Y: 128 Zoom 9: X: 404, Y: 257 Zoom 10: X: 808, Y: 515 Zoom 11: X: 1617, Y: 1031 Zoom 12: X: 3234, Y: 2064

My result is

(TMS false with +1 at Image.WriteTile {tilesY +1}) Zoom 8, X: 404, Y:127 Zoom 9, X: 808, Y: 253 Zoom 10: X: 1617, Y: 505 Zoom 11: X: 3234, Y: 1009 Zoom 12: X: 6468, Y: 2017

(TMS true, Image.WriteTile {tilesY} ) Zoom 8, X: 404, Y:129 Zoom 9, X: 808, Y: 259 Zoom 10: X: 1617, Y: 519 Zoom 11: X: 3234, Y: 1039 Zoom 12: X: 6468, Y: 2079

I can see that it has some geometric progression happening in the Y naming, while X naming is not in correct zoom level.

Gigas002 commented 5 years ago

Thank you for your tests and .tif data.

That MapTiler website uses mercator projection to generate tiles. GTiff2Tiles now is able to create only geodetic (EPSG:4326) tiles.

I've tested Tokyo dataset (approximate location on Google Maps), here, take a look at the results for 11-8 zooms.

  1. MapTiler website, Google tiles. The result should be the same, as MapTiler Pro generates with -mercator and without -tms arguments (see number 3):
Tiles (x_y) Zoom
1819_809 11
909_403 10
454_201 9
227_100 8
  1. MapTiler website, TMS tiles. The result should be the same, as MapTiler Pro generates with -mercator and -tms arguments (see number 4):
Tiles (x_y) Zoom
1819_1241 11
909_620 10
454_310 9
227_155 8
  1. MapTiler Pro with -mercator and without -tms arguments. The result should be the same, as on MapTiler website, Google tiles (see number 1):
Tiles (x_y) Zoom
1819_806 11
909_403 10
454_201 9
227_100 8
  1. MapTiler Pro with -mercator and with -tms arguments. The result should be the same, as on MapTiler website, TMS tiles (see number 2):
Tiles (x_y) Zoom
1819_1241 11
909_620 10
454_310 9
227_155 8

That was all for mercator tiles. These tiles uses a bit different naming and zooming systems, so the values between it and geodesic tiles are differs. Now let's take a look at the geodesic tiles tests results:

  1. Tile.gbif (OpenLayers, uses geodesic non-tms tiles). The results should be same, as GTiff2Tiles generates with tmsCompatible = false and tileY+1 in path (see number 8):
Tiles (x_y) Zoom
3638_618, 3638_619, 3639_618, 3639_619 11
1819_309, 1819_310 10
909_155 9
454_78 8
  1. MapTiler Pro with -geodesic and with -tms arguments. The results should be the same, as GTiff2Tiles generates with tmsCompatible = true and tileY in path (see number 9):
Tiles (x_y) Zoom
3638_1429, 3638_1430, 3639_1429, 3639_1430 11
1819_714, 1819_715 10
909_357 9
454_178 8
  1. MapTiler Pro with -geodesic and without -tms arguments. The results should be the same, as GTiff2Tiles generates with tmsCompatible = false and tileY in path (see number 10):
Tiles (x_y) Zoom
3638_617, 3638_618, 3639_617, 3639_618 11
1819_308, 1819_309 10
909_154 9
454_77 8
  1. GTiff2Tiles with tmsCompatible = false and tileY+1 in path. The results should be the same, as on Tile.gbif (see number 5):
Tiles (x_y) Zoom
3638_618, 3638_619, 3639_618, 3639_619 11
1819_309, 1819_310 10
909_155 9
454_78 8
  1. GTiff2Tiles with tmsCompatible = true and tileY in path. The results should be the same, as MapTiler Pro generates with -geodesic and with -tms arguments (see number 6):
Tiles (x_y) Zoom
3638_1429, 3638_1430, 3639_1429, 3639_1430 11
1819_714, 1819_715 10
909_357 9
454_178 8
  1. GTiff2Tiles with tmsCompatible = false and tileY in path. The results should be the same, as MapTiler Pro generates with -geodesic and without -tms arguments (see number 7):
Tiles (x_y) Zoom
3638_617, 3638_618, 3639_617, 3639_618 11
1819_308, 1819_309 10
909_154 9
454_77 8

As you can see, all the test are passed. I'll test your data too and post the results later, when I'll be able to do so.

By the way, I've tested all this with join algorithm too, and it works the same as crop, with only exception for non-tms tileY+1, which generates some odd output.

Gigas002 commented 5 years ago

I've tested your GeoTiff. There's really have been an issue with empty tiles in output. It occured, because input GeoTiff is 32 bit, and before creating tiles it should be converted to 8 bit (I've added this code with https://github.com/Gigas002/GTiff2Tiles/commit/f08f690f5d08cd604dc0ffa46171fd98d0c4a8ee).

Regarding the tile numbers, it seems to work correctly (I've checked the output data using scheme I've written above, and got the same location on both maptiler and tile.gbif websites).

Feel free to reopen issue and ask more question, if you problem hasn't been solved.

carnegiea commented 5 years ago

after your testing, how is your output in terms of pngsave? In my test it render nothing, unless i changed to jpegsave

Gigas002 commented 5 years ago

Output data seems correct. DId you check it using GTiff2Tiles with the latest commits? Output tiles were transparent because of the issue, that was fixed with today's commit https://github.com/Gigas002/GTiff2Tiles/commit/f08f690f5d08cd604dc0ffa46171fd98d0c4a8ee. On the data, you've send above it now works as it should. Just check out lower levels, you can definetely see the data on 11 zoom for example: 1039