Subscribe to this thread
Home - General / All posts - How to clip an image from a drawing
Corentin
159 post(s)
#16-Jul-20 18:26

Hi,

I am working on the nightlight of Peru and I downloaded a very large image (1.9Go) from NOAA website :

https://ngdc.noaa.gov/eog/viirs/download_dnb_composites.html

The image I downloaded covers half of the southern hemisphere I want to cut (clip) this image so I only have the pixel of Peru.

(Of course I have a drawing with the borders of Peru)

In M8 I would have made a map, selected the shape in the drawing, transfer selection, invert selection and delete the pixels.

How can I do such a thing in M9?

Dimitri


7,413 post(s)
#17-Jul-20 08:20

9 doesn't (yet) have a "clip this raster to the extents of this vector drawing..." command, but I'm sure there's a way to do it now.

One possibility to explore: I haven't actually tried this (if I have time I'll try later today...), but as a thought experiment the following workflow may succeed:

1. Using the methods in this topic, create a new image of the desired number of channels (I'd use RGBA) and resolution, and use Edit - Join to transfer from the Peru drawing into that image opaque values for the alpha channel which the Peru area covers, and transparent values where it does not.

2. You now have an image that covers Peru with opacity where Peru exists and transparency where it does not. Use Edit - Join to successively join the three channels of your very large image into that image. Done.

The above thought experiment is two step: create a starter image that covers Peru at the desired resolution, which is blank except for the alpha channel, and then populate pixel values by joining in from the big image.

As to "desired resolution," I'd probably start by using whatever is the resolution, in terms of the size of a pixel, used by the big image.

Corentin
159 post(s)
#17-Jul-20 11:22

If you confirm that M9 doesn't yet has this functionnality I will stick to my old method for large images of this kind : do it with QGIS (really straightforward there!)

M8 doesn't succeed in loading an image of this size...

Thank you for your answer.

tjhb
10,094 post(s)
#17-Jul-20 22:42

I can help with this. It is not much more complicated in 9 than in 8.

Do you especially need to clip to the irregular shape of Peru, or simply to its bounding box?

Corentin
159 post(s)
#17-Jul-20 23:47

It can be the bounding box...

But I already did it in QGIS. I posted here my problem to know if there was a way of doing it that I missed in the documentation.

If you say it is as easy as in M8 I will try :)

tjhb
10,094 post(s)
#18-Jul-20 01:01

I didn't say that and still don't. I said it is not much more complicated. I hope you will try anyway.

LandSystems73 post(s)
#30-Nov-20 06:12

This is something that I would like to know how to do, ideally using a selection of tiles or the tiles touching a drawing (I appreciated that we can only do full tiles for now).

I have used your merge image technique with the modified rect property (and the re-project method before that), but something like a selection would be really convenient.

I had a go using the join dialogue and a few other options, but I couldn't quite grasp how to set the rect of the output to to the extent of the output image rather than to the extent of the source image with the blank pixels set to be invisible.

rk
621 post(s)
#03-Jun-21 15:07

I need to crop an image today.

~20GB geotiff

Size 122070, 146024 RGBA

COMPRESSION=LZW

38GB manifold project.

I'll try it out.

Attachments:
crop_geotiff2.png

Mike Pelletier

2,122 post(s)
#03-Jun-21 17:50

In case it is helpful, you can manually delete tiles (select in image and then delete in table) not touching the boundary of the desired clip to reduce data size before doing the clip. Also, look in the list of Dimitri's videos for one on clipping images.

rk
621 post(s)
#03-Jun-21 18:15

1. I selected and deleted tiles that had no visible pixels (only [0,0,0,0] pixels). That was more than half of tiles (from 68000+ to 27000+ records in table). Project size reduced less than 2%.

DELETE

FROM 

 [orthomosaic Tiles] 

WHERE

  TileValueMax(TileChannel([Tile], 3)) = 0;

I'm still thinking, how to select tiles that do not touch the area at all. I think I need an "index drawing".

tjhb
10,094 post(s)
#03-Jun-21 18:46

Riivo,

I think this could help as a start.

getting bounding box of tile as a geometry

Tim

rk
621 post(s)
#03-Jun-21 19:05

Magic!

Thank you, Tim.

Riivo

tjhb
10,094 post(s)
#03-Jun-21 19:10

Oops, there were revised versions in this post (same thread).

rk
621 post(s)
#04-Jun-21 14:18

I did it a bit differently.

1. Linked 20GB tiff

2. Dropped it into project (import)

3. Saved project (37GB)

4. I added [Geom] field to [Image Table]. + drawing [Image TileBoxes]

5. Used modified Tim's query (attached) to fill [Geom] field.

6. Spatial Selection in Select pane between [Areas] drawing and [Image TileBoxes] drawing.

7. Reverted selection

8. Tile transform Compose [0, 0, 0, 0] on selection

9. Saved project (13GB)

Exporting tiff probably has a bug. When .ige file reaches 66.7 GB (71,674,637,449 bytes exactly) Manifold says 'Cannot open file.' - will report.

The result is ok, but blocky edges. I will not go into splitting tiles into pixels and back. Or should I try?

Attachments:
blocky_edges.png
TileBoxes.txt

rk
621 post(s)
#03-Jun-21 21:27

Some interesting observations:

Original geotiff opened instantly when linked, thanks to "embedded overviews".

gdalinfo:

Band 1 Block=512x512 Type=Byte, ColorInterp=Red

Overviews: 61035x73012, 30518x36506, 15259x18253, 7630x9127, 3815x4564, 1908x2282, 954x1141, 477x571, 239x286

Mask Flags: PER_DATASET ALPHA

Overviews of mask band: 61035x73012, 30518x36506, 15259x18253, 7630x9127, 3815x4564, 1908x2282, 954x1141, 477x571, 239x286

etc

I deleted all tiles that did not touch the area. About 8000 remained. mapfile size 13GB.

I exported tiff. It took ? time and kind-of-failed with 'Can't read data'. tiff was 6GB, .IGE was 69.99GB. because of the failure, I cannot find export time from log.

The resulting 6gb tiff (COMPRESSION=DEFLATE) , when linked, took ~500 sec to build cache (3gb) before becoming visible.

gdalinfo:

Driver: GTiff/GeoTIFF

Warning 3: VSIFReadL(0000016FFA357020,6*4) @ 4669780 failed in HFAEntry().

No such file or directory

Warning 1: Unable to open orthomosaic-512.aux, it has zero usable bands.

Warning 4: `orthomosaic-512.aux' not recognized as a supported file format.

Warning 3: VSIFReadL(0000016FFA356EC0,6*4) @ 4669780 failed in HFAEntry().

No such file or directory

Warning 1: Unable to open orthomosaic-512.aux, it has zero usable bands.

Warning 4: `orthomosaic-512.aux' not recognized as a supported file format.

Files: orthomosaic-512.tiff

etc 

tiff itself seems OK.

Attachments:
file_sizes.png

tjhb
10,094 post(s)
#03-Jun-21 21:54

Weird? It almost looks as if Manifold is deleting the native-resolution tiles when you instruct it to, but retaining their stored overviews (now disembodied) in its .map database--and then attempting to export them.

Resulting in a large, corrupt TIFF.

Overviews--in any format--should obviously be deleted whenever their parent tiles are deleted.

tjhb
10,094 post(s)
#03-Jun-21 22:16

(That's exactly what you were thinking too Riivo, I expect.)

rk
621 post(s)
#05-Jun-21 20:02

The instant visibility in Manifold of the original tiff was thanks to overviews and tiling (Block=512x512)

After export from Manifold, I used gdal to tile and build overviews.

gdal_translate -of GTiff -co BIGTIFF=YES -co COMPRESS=LZW -co TILED=YES -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 input.tiff output.tiff

~1h

gdaladdo -r lanczos output.tiff

~1h

It is also possible to clip with gdal. good for scripting but much slower than Manifold.

gdalwarp -of GTiff -cutline clip_area.shp -cl clip_area -co BIGTIFF=YES -co COMPRESS=LZW -co TILED=YES -co BLOCKXSIZE=512 -co BLOCKYSIZE=512 input.tiff output.tiff

~18h

gdaladdo -r lanczos output.tiff 

~1h

geozap
264 post(s)
#06-Jun-21 07:11

Images can be clipped by using Georegistration and specifically the Triangulation method.

1. Import the image and set the Projection

2. Create a drawing and set the same Projection. If you already have a drawing with the area you want to keep, use that.

3. On the drawing place Control Points for registration, covering the area from the image you want to clip.

4. Save control points.

5. Open to the image and load the previously saved control points.

6. Open the drawing and on the registration panel select the image to register as as source

7. Reload the control points to the drawing

8. Select triangulation as method

9. Perform registration of the image

You now have a new clipped copy of the image.

geozap
264 post(s)
#06-Jun-21 08:02

An addition:

If the shape that defines the area you want to clip has too many points (as country borders are), it will be probably impractical to manually place possibly hundreds of control points as I wrote in step 3. In that case:

1. Place one control point

2. Save it to create a new Control points drawing

3. Open the drawing with the country limit

4. Run the Convert transformation (convert to point option) to the country drawing

5. Split the created multipoint to points (Split/coordinates)

6. Copy-paste the points at the Control points drawing

7. Create names for the points in the Name field using the Transformation panel (for example you can copy the mdf_id to the Name field using an expression (you have to cast it to nvchar))

You now have the control points ready to be loaded as described in the previous post.

Manifold User Community Use Agreement Copyright (C) 2007-2021 Manifold Software Limited. All rights reserved.