[gdal-dev] gdalbuildvrt Max Value / Pixel Function
Hello, I am looking for functionality to create a VRT of multiple rasters where the VRT pixel value is the maximum of underlying rasters pixel values. I found out that we have pixel functions https://lists.osgeo.org/pipermail/gdal-dev/2016-September/045134.html that we can add to VRTs to get this behavior but just checking to see if there is a direct way to specify this in gdalbuildvrt command line call or do one has to manually inject this in VRT xml. On the same note does the new gdaltileindex provide this functionality of selecting the max value for each pixel out of the box? Regards, Abdul Raheem Siddiqui ERT Inc. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Debugging mask copy issues when using Band interleaved COGs
Owen, Please create a ticket about that. Even Le 13/02/2025 à 04:40, Owen Glofcheski via gdal-dev a écrit : Hey all, I work with Kurt Schwehr, and we've been debugging failures where an explicit mask isn't preserved through a CreateCopy operation to create a COG through the GTiff driver at HEAD. This failure only started occurring after https://github.com/OSGeo/gdal/commit/3bdc81a205bad8342688873e4ed5716f1208e8f5. After the CreateCopy, when we read back the mask, the mask is expanded from the 1-bit mask to the 8-bit mask inside of gtiffoddbitsband.h. However, the mask is not a 1-bit mask when I add logging here; instead, it's the original 8-bit mask (with 255s and 0s). So when it gets expanded to 8-bits, it's no longer the original mask: E.g., 255 0 255 255 becomes 255 255 ... 255 255 0 0 ... 0 0 255 255 ... So my guess is that somewhere along the line we either copy the mask without turning it back into a bitmask, or we lose metadata tracking that the mask is encoded a certain way. Does anything obvious jump out to folks here? The issue *does not* occur when either COPY_SRC_OVERVIEWS="YES" is removed or INTERLEAVE is set to PIXEL (it only fails on band interleaving). Here's a gist containing an example test that fails at head but previously succeeded: https://gist.github.com/orglofch/40cf4d9136978b0384141cc74d1039b6. Thanks, Owen ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev -- http://www.spatialys.com My software is free, but my time generally not. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Problem with warping into existing TIF
Tim, Reviewing the warping logic, I see we have support for a *input* mask, but not creating/updating an *output* mask. AFAICS, there isn't any particular strong reason not to have support for that. Likely just a few hours of coding away. Your best workaround is to generate an output with an alpha band (possibly as a VRT) and use gdal.Translate(dst, src, options="-b 1 -b 2 -b 3 -mask 4") for example to transform a RGBA into a RGB+mask product Even Le 12/02/2025 à 23:47, Tim Harris via gdal-dev a écrit : Not sure if this is a bug, or expected behavior, or user error. If I warp one TIF into another, and both have nodata masks, it seems that gdalwarp isn't updating the destination TIF's nodata mask to unmask the new pixels. Here's a short script to generate two example TIFs. One has a red box in its upper left corner, the other has a green box in its upper right corner. Both have nodata masks so that only the colored boxes are visible. --- import numpy as np from osgeo import gdal gdal.UseExceptions() gdal.SetConfigOption("GDAL_TIFF_INTERNAL_MASK", "YES") drv = gdal.GetDriverByName("GTiff") zeros = np.zeros((3, 1024, 1024)) red = np.copy(zeros) red[0, 100:200, 100:200] = 255 green = np.copy(zeros) green[1, 100:200, 800:900] = 255 ds = drv.Create("red.tif", 1024, 1024, 3, gdal.GDT_Byte) ds.SetGeoTransform((0, 1, 0, 0, 0, -1)) ds.WriteArray(red) ds.CreateMaskBand(gdal.GMF_PER_DATASET) ds.GetRasterBand(1).GetMaskBand().WriteArray(red[0, :, :]) ds = None ds = drv.Create("green.tif", 1024, 1024, 3, gdal.GDT_Byte) ds.SetGeoTransform((0, 1, 0, 0, 0, -1)) ds.WriteArray(green) ds.CreateMaskBand(gdal.GMF_PER_DATASET) ds.GetRasterBand(1).GetMaskBand().WriteArray(green[1, :, :]) ds = None --- Then if you copy red.tif to warp.tif, then warp green.tif into warp.tif: cp red.tif warp.tif gdalwarp green.tif warp.tif I would expect the result to be a combined image with both the red and green boxes. But when I load this result into QGIS, I only see the red box. If I change the layer in QGIS to ignore the mask band, I do see both boxes (and the remaining image is black, as expected). So the pixel data is there, but the mask didn't get updated. Is there a way to get gdalwarp to update the nodata mask automatically? Maybe some hidden setting I'm missing? I can work around this by separately warping the masks, then copying the result into the final TIF's mask band. But it would be nice if the warp took care of this for me. Thanks ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev -- http://www.spatialys.com My software is free, but my time generally not. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] gdalbuildvrt Max Value / Pixel Function
> > Any preliminary documentation on this? > Preliminary, yes: https://gdal--11850.org.readthedocs.build/en/11850/programs/gdal_raster_calc.html (from https://github.com/OSGeo/gdal/pull/11850) Dan ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] gdalbuildvrt Max Value / Pixel Function
WOW! Game changer. Any preliminary documentation on this? On 2/13/25 13:00, Daniel Baston via gdal-dev wrote: the next release of GDAL should have a replacement "gdal_calc" that would take multiple inputs and write a VRT that applies an expression to those inputs ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] gdalbuildvrt Max Value / Pixel Function
Hi Abdul, For the specific case of multiple rasters that share the same extent, the next release of GDAL should have a replacement "gdal_calc" that would take multiple inputs and write a VRT that applies an expression to those inputs ("max", in your case). In the meantime, I think your best option would be to prepare the VRT yourself. Dan On Thu, Feb 13, 2025 at 3:12 PM Abdul Raheem Siddiqui via gdal-dev < gdal-dev@lists.osgeo.org> wrote: > Hello, > > I am looking for functionality to create a VRT of multiple rasters where > the VRT pixel value is the maximum of underlying rasters pixel values. > > I found out that we have pixel functions > https://lists.osgeo.org/pipermail/gdal-dev/2016-September/045134.html > that we can add to VRTs to get this behavior but just checking to see if > there is a direct way to specify this in gdalbuildvrt command line call or > do one has to manually inject this in VRT xml. > > On the same note does the new gdaltileindex provide this functionality of > selecting the max value for each pixel out of the box? > > Regards, > Abdul Raheem Siddiqui > ERT Inc. > > ___ > gdal-dev mailing list > gdal-dev@lists.osgeo.org > https://lists.osgeo.org/mailman/listinfo/gdal-dev > ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] gdalbuildvrt Max Value / Pixel Function
Thanks for your reply, my rasters are not sharing extent but some do overlap. How about GTI, Does that provide this functionality? Do you know. On Thu, Feb 13, 2025, 4:00 PM Daniel Baston wrote: > Hi Abdul, > > For the specific case of multiple rasters that share the same extent, the > next release of GDAL should have a replacement "gdal_calc" that would take > multiple inputs and write a VRT that applies an expression to those inputs > ("max", in your case). In the meantime, I think your best option would be > to prepare the VRT yourself. > > Dan > > On Thu, Feb 13, 2025 at 3:12 PM Abdul Raheem Siddiqui via gdal-dev < > gdal-dev@lists.osgeo.org> wrote: > >> Hello, >> >> I am looking for functionality to create a VRT of multiple rasters where >> the VRT pixel value is the maximum of underlying rasters pixel values. >> >> I found out that we have pixel functions >> https://lists.osgeo.org/pipermail/gdal-dev/2016-September/045134.html >> that we can add to VRTs to get this behavior but just checking to see if >> there is a direct way to specify this in gdalbuildvrt command line call or >> do one has to manually inject this in VRT xml. >> >> On the same note does the new gdaltileindex provide this functionality of >> selecting the max value for each pixel out of the box? >> >> Regards, >> Abdul Raheem Siddiqui >> ERT Inc. >> >> ___ >> gdal-dev mailing list >> gdal-dev@lists.osgeo.org >> https://lists.osgeo.org/mailman/listinfo/gdal-dev >> > ___ gdal-dev mailing list gdal-dev@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/gdal-dev