jblindsay / whitebox-tools

An advanced geospatial data analysis platform
https://www.whiteboxgeo.com/
MIT License
949 stars 161 forks source link

clip_raster_to_polygon issue #102

Closed mlavy closed 4 years ago

mlavy commented 4 years ago

the clip_raster_to_polygon tool fail to execute returning:

**********************************
* Welcome to ClipRasterToPolygon *
**********************************
Reading data...
thread 'main' panicked at 'capacity overflow', src\liballoc\raw_vec.rs:777:5
note: run with `RUST_BACKTRACE=1` environment variable to display a backtrace.

Process finished with exit code 0

the tool was used as following:

wbt.clip_raster_to_polygon(
    i=input_raster,
    polygons=task_boundary,
    output=out_clipped_raster,
    maintain_dimensions=False,
    #callback=default_callback
)

attached the input data in Geographic CRS: 4326

issue.zip

mlavy commented 4 years ago

everything is attached in the github issue. Anyway here attached again

On Thu, Jul 23, 2020 at 3:59 PM John Lindsay notifications@github.com wrote:

I'll need more information to go on than that to diagnose this problem. Can you share your raster and polygon files with me?

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/jblindsay/whitebox-tools/issues/102#issuecomment-663023518, or unsubscribe https://github.com/notifications/unsubscribe-auth/AMWMWHTGY6L6FKZULBZP5X3R5A64VANCNFSM4PFX6HRQ .

--

Muriel Lavy

e-mail: murielavy@g muriel.lavy@polito.itmail.com

Skype: muriel.lavy

jblindsay commented 4 years ago

Yeah, I'm sorry I didn't see that until after I posted my comment, which I then deleted. I noticed that it runs if you check 'Maintain input raster dimensions' but the output raster is empty. I'll have to dig a little deeper into the issue. It likely has something to do with the use of geographic coordinates instead of a projected CRS.

mlavy commented 4 years ago

It seems to work with the "maintain_dimensions=True" but the result is an empty raster with a wrong spatial extension. Thanks for checking

On Thu, Jul 23, 2020 at 4:04 PM John Lindsay notifications@github.com wrote:

Yeah, I'm sorry I didn't see that until after I posted my comment, which I then deleted. I noticed that it works if you check 'Maintain input raster dimensions'. I'll have to dig a little deeper into the issue.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/jblindsay/whitebox-tools/issues/102#issuecomment-663026428, or unsubscribe https://github.com/notifications/unsubscribe-auth/AMWMWHU26OSPFZR27IGDU4TR5A7QRANCNFSM4PFX6HRQ .

--

Muriel Lavy

e-mail: murielavy@g muriel.lavy@polito.itmail.com

Skype: muriel.lavy

jblindsay commented 4 years ago

So I think it comes down to your GeoTIFF--I'm not sure if I've told you before just how much I dislike the GeoTIFF format but it is truly the bane of my existence. I spend an inordinate amount of my time trying to resolve issues with the GeoTIFF reader/writer in WBT. The GeoTIFF file format is all together too flexible and complex. I would ban it from geoscience if I could.

Anyhow, mandatory GeoTIFF rant aside, running your input file through the PrintGeoTiffTags tool shows the following:

*******************************
* Welcome to PrintGeoTiffTags *
*******************************
Little-endian byte order used
Classic TIFF format
TIFF Tags:
Number of tags: 19
ImageWidth (code=256 type=DT_Short): 919
ImageLength (code=257 type=DT_Short): 728
BitsPerSample (code=258 type=DT_Short): 32
Compression (code=259 type=DT_Short): "Deflate" (8)
PhotometricInterpretation (code=262 type=DT_Short): "BlackIsZero" (1)
SamplesPerPixel (code=277 type=DT_Short): 1
PlanarConfiguration (code=284 type=DT_Short): "Contiguous" (1)
Predictor (code=317 type=DT_Short): "None" (1)
TileWidth (code=322 type=DT_Short): 256
TileLength (code=323 type=DT_Short): 256
TileOffsets (code=324 type=DT_Long count=12 offset=290): [504, 62533, 123578, 182615, 218269, 283388, 347269, 411217, 450160, 499571, 547277, 596099]
TileByteCounts (code=325 type=DT_Long count=12 offset=242): [62029, 61045, 59037, 35654, 65119, 63881, 63948, 38943, 49411, 47706, 48822, 30051]
SampleFormat (code=339 type=DT_Short): "Floating point data" (3)
ModelPixelScaleTag (code=33550 type=DT_Double count=3 offset=344): [0.00001035316313322134, 0.00001035316313322134, 0.0]
ModelTiepointTag (code=33922 type=DT_Double count=6 offset=368): [0.0, 0.0, 0.0, -101.3756775748901, 38.39447770745262, 0.0]
GeoKeyDirectoryTag (code=34735 type=DT_Short count=32 offset=416): [1, 1, 0, 7, 1024, 0, 1, 2, 1025, 0, 1, 1, 2048, 0, 1, 4326, 2049, 34737, 7, 0, 2054, 0, 1, 9102, 2057, 34736, 1, 1, 2059, 34736, 1, 0]
GeoDoubleParamsTag (code=34736 type=DT_Double count=2 offset=480): [298.257223563, 6378137.0]
GeoAsciiParamsTag (code=34737 type=DT_ASCII count=8 offset=496): WGS 84|
GDAL_NODATA (code=42113 type=DT_ASCII count=6 offset=338): -9999
Number of tags: 15
NewSubFileType (code=254 type=DT_Long offset=1): 1
ImageWidth (code=256 type=DT_Short): 460
ImageLength (code=257 type=DT_Short): 364
BitsPerSample (code=258 type=DT_Short): 32
Compression (code=259 type=DT_Short): "Deflate" (8)
PhotometricInterpretation (code=262 type=DT_Short): "BlackIsZero" (1)
SamplesPerPixel (code=277 type=DT_Short): 1
PlanarConfiguration (code=284 type=DT_Short): "Contiguous" (1)
Predictor (code=317 type=DT_Short): "None" (1)
TileWidth (code=322 type=DT_Short): 128
TileLength (code=323 type=DT_Short): 128
TileOffsets (code=324 type=DT_Long count=12 offset=626384): [627046, 674173, 720549, 764485, 791171, 840541, 888363, 936341, 965527, 1002583, 1037337, 1073762]
TileByteCounts (code=325 type=DT_Long count=12 offset=626336): [47127, 46376, 43936, 26686, 49370, 47822, 47978, 29186, 37056, 34754, 36425, 22647]
SampleFormat (code=339 type=DT_Short): "Floating point data" (3)
GDAL_NODATA (code=42113 type=DT_ASCII count=6 offset=626432): -9999
Number of tags: 15
NewSubFileType (code=254 type=DT_Long offset=1): 1
ImageWidth (code=256 type=DT_Short): 230
ImageLength (code=257 type=DT_Short): 182
BitsPerSample (code=258 type=DT_Short): 32
Compression (code=259 type=DT_Short): "Deflate" (8)
PhotometricInterpretation (code=262 type=DT_Short): "BlackIsZero" (1)
SamplesPerPixel (code=277 type=DT_Short): 1
PlanarConfiguration (code=284 type=DT_Short): "Contiguous" (1)
Predictor (code=317 type=DT_Short): "None" (1)
TileWidth (code=322 type=DT_Short): 128
TileLength (code=323 type=DT_Short): 128
TileOffsets (code=324 type=DT_Long count=4 offset=626640): [1096409, 1145793, 1184618, 1203508]
TileByteCounts (code=325 type=DT_Long count=4 offset=626624): [49384, 38825, 18890, 15700]
SampleFormat (code=339 type=DT_Short): "Floating point data" (3)
GDAL_NODATA (code=42113 type=DT_ASCII count=6 offset=626656): -9999
Number of tags: 15
NewSubFileType (code=254 type=DT_Long offset=1): 1
ImageWidth (code=256 type=DT_Short): 115
ImageLength (code=257 type=DT_Short): 91
BitsPerSample (code=258 type=DT_Short): 32
Compression (code=259 type=DT_Short): "Deflate" (8)
PhotometricInterpretation (code=262 type=DT_Short): "BlackIsZero" (1)
SamplesPerPixel (code=277 type=DT_Short): 1
PlanarConfiguration (code=284 type=DT_Short): "Contiguous" (1)
Predictor (code=317 type=DT_Short): "None" (1)
TileWidth (code=322 type=DT_Short): 128
TileLength (code=323 type=DT_Short): 128
TileOffsets (code=324 type=DT_Long offset=1219208): 1219208
TileByteCounts (code=325 type=DT_Long offset=31593): 31593
SampleFormat (code=339 type=DT_Short): "Floating point data" (3)
GDAL_NODATA (code=42113 type=DT_ASCII count=6 offset=626848): -9999
Number of tags: 15
NewSubFileType (code=254 type=DT_Long offset=1): 1
ImageWidth (code=256 type=DT_Short): 58
ImageLength (code=257 type=DT_Short): 46
BitsPerSample (code=258 type=DT_Short): 32
Compression (code=259 type=DT_Short): "Deflate" (8)
PhotometricInterpretation (code=262 type=DT_Short): "BlackIsZero" (1)
SamplesPerPixel (code=277 type=DT_Short): 1
PlanarConfiguration (code=284 type=DT_Short): "Contiguous" (1)
Predictor (code=317 type=DT_Short): "None" (1)
TileWidth (code=322 type=DT_Short): 128
TileLength (code=323 type=DT_Short): 128
TileOffsets (code=324 type=DT_Long offset=1250801): 1250801
TileByteCounts (code=325 type=DT_Long offset=8568): 8568
SampleFormat (code=339 type=DT_Short): "Floating point data" (3)
GDAL_NODATA (code=42113 type=DT_ASCII count=6 offset=627040): -9999
GeoKey Info:
Version: 1
Key revision: 1.0
GTModelTypeGeoKey (code=1024): ModelTypeGeographic (2)
GTRasterTypeGeoKey (code=1025): RasterPixelIsArea (1)
GeographicTypeGeoKey (code=2048): GCS_WGS_84 (GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433],AUTHORITY["EPSG",4326]])
GeogCitationGeoKey (code=2049, type=ASCII): WGS 84
GeogAngularUnitsGeoKey (code=2054): Angular_Degree (9102)
GeogSemiMajorAxisGeoKey (code=2057, type=Double, count=1): 6378137
GeogInvFlatteningGeoKey (code=2059, type=Double, count=1): 298.257223563

I'm not sure how you created this GeoTIFF, but I've never seen one like it before now. There are several 'sub images' (I suspect that it has had pyramiding applied, but I'd be curious to know if this is the case), the last one having rows x columns of 58x46. Therein lies the problem. WBT is seeing the images number of rows and columns as being this last sub-image and when it calculates the rows/columns of positions within the polygon, all of them fall outside of the predicted raster location. It's a problem.

I'll see if I can make another patch to WBT's GeoTIFF reader to only read the first of multiple sub-images. In the meantime, iIf you convert your file to a BIL, the tool works as expected.

jblindsay commented 4 years ago

I've just committed code to fix the GeoTIFF reader to only read the first of multiple sub images contained in a GeoTIFF. This resolves the issue that you ran into. Unfortunately, I just released v1.3.1 this morning and so you will either need to wait until the next release for this new feature, or will need to compile from source.

mlavy commented 4 years ago

The GeoTIFF was generated using GDAL (gdal_translate) with pyramids containing TILES of 256 x 256 ( block size) and compression: DEFLATE of 9 levels. Inspecting the raster in QGIS you can actually see that the last sub imaged is 58x46.

image