Opened 12 years ago

Closed 12 years ago

Last modified 12 years ago

#1839 closed defect (fixed)

[raster] raster2pgsql is unable to handle TIFF subdatasets

Reported by: turbothomas Owned by: Bborie Park
Priority: medium Milestone: PostGIS 2.0.1
Component: raster Version: 2.0.x
Keywords: raster2pgsql, TIFF, subdataset, history Cc: rouault

Description

When loading TIFFs with subdatasets, only the first subdataset is loaded, all others are ignored.

Attachments (1)

raster2pgsqlx64_r9822.zip (442.3 KB ) - added by robe 12 years ago.
2.0.1SVN raster2pgsql 64-bit PostgreSQL windows

Download all attachments as: .zip

Change History (29)

comment:1 by Bborie Park, 12 years ago

Component: loader/dumperraster
Owner: changed from mcayland to Bborie Park
Summary: raster2pgsql is unable to handle TIFF subdatasets[raster] raster2pgsql is unable to handle TIFF subdatasets

comment:2 by Bborie Park, 12 years ago

Can you provide the output from running gdalinfo on the tiff file?

in reply to:  2 comment:3 by turbothomas, 12 years ago

Replying to dustymugs:

Can you provide the output from running gdalinfo on the tiff file?

Driver: GTiff/GeoTIFF
Files: pathto/imagename.tif
Size is 611, 1
Coordinate System is `'
Metadata:
  TIFFTAG_IMAGEDESCRIPTION=CXmlTiff type 
  TIFFTAG_RESOLUSTIONUNIT=2 (pixels/inch)
  TIFFTAG_XRESOLUTION=100
  TIFFTAG_YRESOLUTION=100
ImageStructure Metadata:
  INTERLEAVE=BAND
Subdatasets:
  SUBDATASET_1_NAME=GTIFF_DIR:1:pathto/imagename.tif
  SUBDATASET_1_DESC=Page 1 (611P x 1L x 1B)
  SUBDATASET_2_NAME=GTIFF_DIR:2:pathto/imagename.tif
  SUBDATASET_2_DESC=Page 1 (392P x 514352L x 1B)
Corner Coordinaes:
Upper Left  (   0.0, 0.0)
Lower Left  (   0.0, 1.0)
Upper Right ( 611.0, 0.0)
Lower Right ( 611.0, 1.0)
Center      ( 305.5, 0.5)
Band 1 Block=611x1 Type=Byte, ColorInterp=Gray

comment:4 by Bborie Park, 12 years ago

Status: newassigned

Should be fixed in r9807. Please close this ticket once you've tested the subdataset. I don't have any GeoTIFF files with subdatasets to test with.

comment:5 by robe, 12 years ago

Bborie,

Do you need to make the same change to trunk or is trunk not have the issue?

turbothomas,

Sounds like you are testing on Windows maybe. Let me know if you need a build to test with. I was going to get one out this week anyway and can probably expedite it a bit or at least give you a new raster2pgsql to test now before I release.

comment:6 by Bborie Park, 12 years ago

Yes, I do need to make the same change to trunk. Just holding off on it until I can resolve #1840.

comment:7 by Bborie Park, 12 years ago

Keywords: history added

Fixed for 2.1 in r9809

comment:8 by rouault, 12 years ago

Cc: rouault added

If I'm allowed to comment a bit about the code, I don't think that the approach taken in isSubDataset() is so great in term of exaustiveness and maintainability. I'm not sure you have captured all the GDAL drivers that support subdatasets (actually I know some others, and the list can potentially grow with newer GDAL version), and it should not be the job of PostGIS to know their syntax. I'd recommand using GDALIdentify() and if it succeeds, then it's a GDAL dataset or subdataset.

Another related remark about the code around line 2537 that uses fopen() to test if it is a raster, that is based on fopen(). GDAL supports opening files on exotic filesystems like /vsizip/ , /vsicurl/ or even things that are not files, like a PostgisRaster dataset, or a GeoRaster dataset… Same as above, GDALIdentify() would be better.

So all in all I think you could just drop isSubDataset() and use GDALIdentify().

(Sorry if I'm perhaps misunderstanding the code. My comments are only based on code reading.)

comment:9 by Bborie Park, 12 years ago

Thanks rouault! I'll definitely look into using GDALIdentify() instead.

in reply to:  5 comment:10 by turbothomas, 12 years ago

Replying to robe:

yes, I am working with Windows 7 64-bit. If I could get a build, I can test on my datasets.

comment:11 by Bborie Park, 12 years ago

turbothomas,

Is any of the TIFF dataset publicly available? I'd like to be able to have an example of one for future testing…

comment:12 by rouault, 12 years ago

try http://svn.osgeo.org/gdal/trunk/autotest/gcore/data/twoimages.tif

or use tiffcp, like in : tiffcp src.tif src.tif dst.tif

comment:13 by turbothomas, 12 years ago

If that doesn't do I can upload one on a file server tomorrow.

comment:14 by robe, 12 years ago

turbothomas,

Are you using the 64-bit PostgreSQL or 32-bit. Sorry been busy this week and will be much of weekend trying to meet our book deadline.

Attached is the 64-bit raster2pgsql from 2.0.1 branch. Should be a drop in replace for what you've got already. Let me know if you need the 32-bit one.

by robe, 12 years ago

Attachment: raster2pgsqlx64_r9822.zip added

2.0.1SVN raster2pgsql 64-bit PostgreSQL windows

in reply to:  14 comment:15 by turbothomas, 12 years ago

I am using 64-bit PostgreSQL. Thanks, I will test it.

in reply to:  14 comment:16 by turbothomas, 12 years ago

Only got to test it now.

The tiling seems to work fine, now. For the subsets, I am still getting the same results. Only the first one is loaded, the others are ignored. Do I have to set a special flag?

I uploaded two files I am using. One with two subsets, the other one with several

http://www.sendspace.com/file/8u3vzq http://www.sendspace.com/file/lzqdog

comment:17 by Bborie Park, 12 years ago

What are you passing to raster2pgsql? If I'm correct based upon with working with other raster formats (HDF4, HDF5, NetCDF) with subdatasets, you can't pass the raster file in as that automatically assumes the first subdataset. Instead, you need to pass in something like:

SUBDATASET_1_NAME=GTIFF_DIR:1:pathto/imagename.tif

The above is from the gdalinfo output you posted in a prior comment. So, the following only gets the first subdataset.

raster2pgsql -s 0 -t 50x50 pathto/imagename.tif tablename

The following gets the specific subdataset

raster2pgsql -s 0 -t 50x50 SUBDATASET_1_NAME=GTIFF_DIR:1:pathto/imagename.tif tablename

The behavior is as such because a raster can have multiple subdatasets, each of which can have multiple bands. Actually, the behavior is exactly the same as with any of the gdal utilities.

in reply to:  17 comment:18 by turbothomas, 12 years ago

ok, I didn't get that part with the name and subdataset.

Tested it now, passing the name as dustymugs described above.

Got back

ERROR: Unable to read raster file: SUBDATASET_2_NAME=GTIFF_DIR:2:myfile.tif

(same for the first subdataset)

comment:19 by Bborie Park, 12 years ago

Can you run the following?

gdalinfo SUBDATASET_2_NAME=GTIFF_DIR:2:myfile.tif

If that doesn't work, raster2pgsql won't work either.

Also, always post your raster2psql string as there is no way for anyone helping to replicate your situation.

in reply to:  19 comment:20 by turbothomas, 12 years ago

sorry, before I was trying

raster2pgsql SUBDATASET_2_NAME=GTIFF_DIR:2:myfile.tif  tablename

Now I tried

gdalinfo SUBDATASET_2_NAME=GTIFF_DIR:2:myfile.tif 

and got

ERROR 4: 'SUBDATASET_2_NAME=GTIFF_DIR:2:myfile.tif ' does not exist in the file system, and is not recognised as a supported dataset name.

gdalinfo failed - unable to open 'SUBDATASET_2_NAME=GTIFF_DIR:2:myfile.tif '.

comment:21 by Bborie Park, 12 years ago

Now that I look at it, you're not supposed to use the "SUBDATASET_2_NAME=" part. Only use the "GTIFF_DIR:2:myfile.tif" after the equals sign. I see that one of my prior comments let you astray. Sorry about that. So…

raster2pgsql GTIFF_DIR:2:myfile.tif  tablename

in reply to:  21 comment:22 by turbothomas, 12 years ago

ok,

gdalinfo GTIFF_DIR:2:myfile.tif

Driver GTIFF/GeoTIFF
Files: myfile.tif
Size is 349, 77520
Coordinate System is `´
Metadata:
  TIFFTAG_DOCUMENTNAME=MICROTEC3DFAST
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  TIFFTAG_XRESOLUTION=100
  TIFFTAG_YRESOLUTION=100
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left (   0.0,    0.0)
Lower Left (   0.0,77520.0)
Upper Right( 340.0.    0.0)
Lower Right( 340.0,77520.0)
Center     ( 170.0,38760.0)
Band 1 Block=340x8 Type=Int16, ColorInterp=Gray

That works :)

Than

raster2pgsql GTIFF:2:myfile.tif test > test.sql

Processing 1/1: GTIFF_DIR:2:myfile.tif
INFO: Using default geotransform matrix (0, 1, 0, 0, 0, -1) for raster: GTIFF_DIR:2:myfile.tif
ERROR: rt_band_set_pixel_line: Coordinates out of range
ERROR: rt_band_set_pixel_line: Coordinates out of range
ERROR: rt_band_set_pixel_line: Coordinates out of range
...

That doesn't work :(

comment:23 by Bborie Park, 12 years ago

Now that is strange. I'll have to do some testing with the same rasters you've provided.

comment:24 by Bborie Park, 12 years ago

Well, that was fun. PostGIS Raster has a maximum permitted width/height size of 65535 x 65535. Unfortunately, there was no obvious error message when you attempt to exceed that limit. I was testing your Eiche*.tif file and hit it using…

raster2pgsql GTIFF_DIR:2:Eiche_1_1_2011_08_17_768_M2_3D_5mm.tif blabla 

So, the lesson learned is that I need to add an explicit error message when someone hits this maximum and you really shouldn't call it the way I did above. Please tile your raster to be at most 65535 x 65535.

And before anyone asks why the maximum width and height is set to 65535 x 65535, imagine the size of the following…

65535 width * 65535 height * 8 bytes per pixel (64BF) = 34358689800 bytes (not counting overhead)

That would be a very monstrous number considering the maximum size permitted in a single field of PostgreSQL is 1 GB (http://www.postgresql.org/about/) .

in reply to:  24 comment:25 by turbothomas, 12 years ago

Well, for the testing I was just to lazy to type in the tiling flag. Normally I would have used it.

But it seems that we could identify another small bug that way.

Thanks for the help

Thomas

comment:26 by Bborie Park, 12 years ago

Resolution: fixed
Status: assignedclosed

Added appropriate error message for 2.0 in r9839. For trunk (2.1) in r 9838.

comment:27 by pracine, 12 years ago

I would have think the message to be in raster2pgsql also. No?

comment:28 by Bborie Park, 12 years ago

Nope. The error occurs in the C API, specifically in rt_raster_new().

Note: See TracTickets for help on using tickets.