[gdal-dev] Help with band pixel functions

[Adriana Parra]

Hello Steven,

Thank you for your reply. I should have mentioned that my ultimate goal is
to mosaic more than two images, so I don't know if the workflow that you
suggest would work in that case. I did use your suggestion of masking
values, and changed the VRT python function to this:

import numpy as np

def average(in_ar, out_ar, xoff, yoff, xsize, ysize,
raster_xsize,raster_ysize, buf_radius, gt, **kwargs):
    x = np.ma.masked_greater(in_ar, 100)
    np.mean(x, axis = 0,out = out_ar, dtype = 'uint8')


That worked, but now I have a new issue. I want to try a maximum function
instead of mean to evaluate if the resulting mosaic would look seamless, I
tried this:

import numpy as np

def average(in_ar, out_ar, xoff, yoff, xsize, ysize,
raster_xsize,raster_ysize, buf_radius, gt, **kwargs):
    x = np.ma.masked_greater(in_ar, 100)
    out_ar[:] = np.ma.max(x, axis = 0, fill_value=0)

The max function works fine in areas where there is data, but for some
reason, areas of overlapping NA values result in a pixel value of 63. I
tried setting the fill value of the masked areas to zero but that didn't
work. I appreciate any suggestions.

Kind regards,

Adriana.

On Mon, Jun 24, 2019 at 7:38 PM Adriana Parra <adrianaparra at nevada.unr.edu>
wrote:

> Hello,
>
> I have been trying to mosaic two images that partially overlap using band
> pixel functions in Python. My images contain unsigned integers (uint8) and
> the no data value is 255. I considered using the nanmean function so that I
> would get average values in overlapping areas of the two images and
> maintain the original value of the non-overlapping areas, but for some
> reason the NAN values are being interpreted as zero, resulting in half of
> the original value in non-overlapping areas.  I would really appreciate any
> feedback. I am using GDAL 3.1.0 and python 3.6. Bellow is the code I used
> and I am sending attached the input files.
> gdalbuildvrt  mosaic1.vrt -input_file_list  input_files.txt -srcnodata
> "255"
>
>  import numpy as np
> def average(in_ar, out_ar, xoff, yoff, xsize, ysize,
> raster_xsize,raster_ysize, buf_radius, gt, **kwargs):
>     np.round_(np.clip(np.nanmean(in_ar, axis = 0, dtype = 'uint8'),0,255),
>               out = out_ar)
>
> gdal_translate --config GDAL_VRT_ENABLE_PYTHON YES mosaic1.vrt mosaic1.tif
>
>  Mosaic1.jpeg
> <https://drive.google.com/a/nevada.unr.edu/file/d/1YgugcVnU6d5BmPYkYcotqBi3fFQcPpY8/view?usp=drive_web>
>  raw_files.jpeg
> <https://drive.google.com/a/nevada.unr.edu/file/d/1SQoH58V3MGtP-oXYnPYsf66snFOsWagC/view?usp=drive_web>
>  test_files.7z
> <https://drive.google.com/a/nevada.unr.edu/file/d/14No5nblxhNKeOU8k86Gfan7utYIUSqRU/view?usp=drive_web>
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20190715/86b646df/attachment.html>