libvips / pyvips

python binding for libvips using cffi
MIT License
636 stars 49 forks source link

Slice a huge tif fails - unable to call VipsForeignSaveTiffFile #487

Closed Hutch07 closed 2 months ago

Hutch07 commented 2 months ago

I am trying to open a very large tif file 15GiB and slice out an area of interest. Then I am trying to write the slice to a new tif file. The following error pops up on the write.

Exception has occurred: Error unable to call VipsForeignSaveTiffFile wbuffer_write: write failed windows error: The printer is out of paper.

This is a Windows 11 machine: Using Python 3.11.

   image_vp = pyvips.Image.new_from_file(image_path)
    image_w = image_vp.width
    image_h = image_vp.height
    sliced_image = image_vp.crop(pix_img_offsetX, pix_img_offsetY, clip_pix_width, clip_pix_height)
     # Save the sliced image
      fname = "clip_" + name_list[i]+ ".tif"
      out_path = str(Path(out_dir, fname))
      print(sliced_image.width, "Width.  Trying to write to ",out_path)
      sliced_image.write_to_file(out_path)`
jcupitt commented 2 months ago

Hi @Hutch07,

The error message is a little misleading heh, you're probably out of disc.

If you run:

   image_vp = pyvips.Image.new_from_file(image_path)

libvips will try to open your image for random access. If your TIFF file does not support random access (plain strip TIFFs do not, for example), then libvips will convert the image to a random access format in a temporary file. This can need a lot of disc space.

The best solution is to use a TIFF variant which supports random access -- interleaved tiled TIFF is best.

If you're unable to use random access TIFF, then I would do:

image = pyvips.Image.new_from_file(image_path, access="sequential")
tile = image.crop(left, top, width, height)
tile.write_to_file("tile.png")

Now libvips will fetch pixels from the file and throw them away until it reaches scanline top, then fetch pixels until the tile is full, then write the tile and stop.

You will only be able to crop out tiles in scanline order, plus it will be fast at the top and very slow at the bottom. But it won't use gigabytes of disc, phew.

You can convert a strip TIFF to a tiled TIFF with eg.:

vips copy huge.tif tiled-huge.tif[tile,compression=jpeg]

It might be a better solution. What are you trying to achieve?

Hutch07 commented 2 months ago

Thanks for the help. I have 275 Regions Of Interest that I want to crop out of a huge geotiff and then create smaller geotiffs to display on a leaflet map. The Geotiff is Width: 833952, Height: 177017. I was going to write a method to PIL.getdata() to get a single line of 833952 and sort them into the 275 ROI's. This might be better. I was originally trying to GDAL tile the Geotiff into TMS/Google tiles but it is taking forever and I only care about a much smaller area.

jcupitt commented 2 months ago

Ah OK. I would just convert directly to a single, huge leaflet map with eg.

vips dzsave huge.tif my-leaflet-tiles --layout google

It'll write all the tiles into my-leaflet-tiles/. It's reasonably quick, and won't use much memory or disc, except for the huge number of tiles heh. Maybe test on a smaller image first though!

jcupitt commented 2 months ago

I tried on a slide file I have:

$ vipsheader f0c003e9-2148-aa33-8135-008f5f76698c_215017.svs
f0c003e9-2148-aa33-8135-008f5f76698c_215017.svs: 183237x40271 uchar, 4 bands, srgb, openslideload
$ /usr/bin/time -f %M:%e vips dzsave f0c003e9-2148-aa33-8135-008f5f76698c_215017.svs[rgb] x --layout google
2686732:157.37

(don't use the rgb flag, that's just for openslide)

It converted in 2.5m and 2.7gb of memory. Your image is about 20x larger, so it should convert in under an hour and need perhaps 10gb of ram (memory use scales with horizontal image size).

Hutch07 commented 2 months ago

I did not realize that this was possible with pyvips. I will start reading info for this and trying it out shortly.

libvips pyramids

jcupitt commented 2 months ago

Yes, you need:

https://www.libvips.org/API/current/VipsForeignSave.html#vips-dzsave

You can use it from python with eg.:

image = pyvips.Image.new_from_file(image_path, access="sequential")
image.dzsave("my-tile-dir", layout="google")
Hutch07 commented 2 months ago

Seems to be making tiles but the levels and google/slippy tile names are off. The geotiff upper left corner is at Latitude = 41.3094838551479 Longitude = -88.45198898587637 XY @EPSG:6455 US Feet = x 951651.556 , y 1691000.711

This looks to be zoom level 22 so I should have up to Directory 22. The zoom level directory should contain subdirectories with the EPSG:3857 X dimension and files with the Y dimension. I think I can convert these but is there something I can do to get them in the correct format before tiling?

Are the upper left coords coming from the Geotiff or a .tfw world file? Or should I be calculating the upper left? upper left tile directory at Zoom Level 22 should be tile 1066611 with file name 2088751.png at zoom 22

jcupitt commented 2 months ago

libvips doesn't know about any of the geotiff tags, so you always get the top-left tile as (0, 0). You'd need to write a little python to rename them all.

Hutch07 commented 2 months ago

Got it. I think I will have to scale the images also since the tiff can be any scale and the google tiles must be the exact meters per pixel for the scale level. Thanks

jcupitt commented 2 months ago

You can scale and save in one operation (I expect you know) eg.:

image = pyvips.Image.new_from_file(image_path, access="sequential")
image = image.resize(1.123564)
image.dzsave("my-tile-dir", layout="google")

Will scale up by 1.123564x during the save.

Hutch07 commented 2 months ago

Does pyvips tile like TMS? Zoom/X/Y.png? I had great success with dzsave but the tiles look like they may be Zoom/Y/X.png.


image_out.dzsave("C:\pyvip_test\Test2", layout="google", background=[0], depth="one") ```
jcupitt commented 2 months ago

No, dzsave is google-maps style, so zoom/y/x.suffix. You'd need to write a little python to move the files around.

Hutch07 commented 2 months ago

I see the layout now. I can just change the leaflet website to use the google style - https://www.microimages.com/documentation/TechGuides/78googleMapsStruc.pdf Thanks again.

Hutch07 commented 2 months ago

Is there a way to export the tiles as .png? Found it

image_out.dzsave("C:\pyvip_test\Test2", layout="google", background=[0], suffix='.png', depth="one")
jcupitt commented 2 months ago

You can pass params to the saver as well, eg. suffix=".png[interlace]".