Updating attributes with band math

[mmann1123]

Maybe a minor tweek to improve band math.

import geowombat as gw
from geowombat.data import rgbn_suba, rgbn_subb

with gw.config.update(sensor='rgbn', scale_factor=0.0001):
    with gw.open([rgbn_suba, rgbn_subb],
          band_names=['blue', 'green', 'red', 'nir'],
          mosaic=True,
          bounds_by='union') as ds_mos:

              # calculate NDVI
              ndvi = ds_mos.gw.ndvi()
              print(ndvi)

returns

<xarray.DataArray (band: 1, y: 283, x: 449)>
dask.array<broadcast_to, shape=(1, 283, 449), dtype=float64, chunksize=(1, 64, 64), chunktype=numpy.ndarray>
Coordinates:
  * y        (y) float64 2.05e+06 2.05e+06 2.05e+06 ... 2.049e+06 2.049e+06
  * x        (x) float64 7.929e+05 7.929e+05 7.929e+05 ... 7.952e+05 7.952e+05
  * band     (band) <U4 'ndvi'
Attributes:
    transform:    (5.0, 0.0, 792925.0, 0.0, -5.0, 2050115.0)
    crs:          +init=epsg:32618
    res:          (5.0, 5.0)
    is_tiled:     0
    nodatavals:   (0.0, 0.0, 0.0, 0.0)
    scales:       (1.0, 1.0, 1.0, 1.0)
    offsets:      (0.0, 0.0, 0.0, 0.0)
    resampling:   nearest
    pre-scaling:  0.0001
    sensor:       rgbn
    drange:       (-1, 1)

Shouldn't attributes be updated to reflect a single band after running band math? For instance nodatavals: (0.0), scales: (1.0) etc? This can be problematic I think when going to write out NDVI and needing to reset attributes.

[mmann1123]

Maybe I should show the exact issue I am seeing, maybe not related to attributes.

with gw.config.update(sensor='rgbn', scale_factor=0.0001):
    with gw.open([rgbn_suba, rgbn_subb],
          band_names=['blue', 'green', 'red', 'nir'],
          mosaic=True,
          bounds_by='union') as ds_mos:
                
                ndvi = ds_mos.gw.ndvi()

                # Xarray drops attributes
                attrs = ndvi.attrs.copy()
                ndvi.attrs['nodatavals'] = (0.0)
                ndvi.attrs['scales'] = (1.0)
                ndvi.attrs['offsets'] = (0.0)
 
                ndvi.attrs = attrs
                print(attrs)
                # Write the data to a GeoTiff
                ndvi.gw.to_raster( '/home/mmann1123/Desktop/output2.tif',
                                    verbose=1,
                                    n_workers=4,    # number of process workers sent to ``concurrent.futures``
                                    n_threads=2,    # number of thread workers sent to ``dask.compute``
                                    n_chunks=200)   # number of window chunks

throws:

ValueError: Source shape (1, 64) is inconsistent with given indexes 1

I am assuming it has to do with the use of the config, but not sure there is a config for single band...

[jgrss]

I agree, the attributes should be updated for band math.

I will see if I can figure out the shape issue. One quick note--if you specify the sensor in the configuration then you don't have to specify the band names in gw.open(). But, I see that I named the files incorrectly. They are stored as blue, green, red, NIR, not red, gree, blue, NIR.

[jgrss]

Ok, the attributes were not the issue when writing the results to file. The problem seems to be the chunk size. The image block size is 64x64, which led to a single row/column at the image edge. Using a larger chunk size fixed it.

>>> with gw.config.update(sensor='bgrn', scale_factor=0.0001):
>>>     
>>>     with gw.open([rgbn_suba, rgbn_subb],
>>>                  mosaic=True,
>>>                  bounds_by='union',
>>>                  chunks=256) as ds_mos:
>>> 
>>>         attrs = ds_mos.attrs.copy()
>>>         
>>>         # Calculate the NDVI
>>>         ndvi = ds_mos.gw.ndvi()
>>>         
>>>         # Mask the NDVI
>>>         masker = ndvi >= -0.25
>>>         ndvi = ndvi * masker
>>>         
>>>         # The 'no data' values are set as None from `ds_mos.gw.ndvi()`
>>>         # After masking, set the 'no data' values as 0.
>>>         ndvi.attrs = attrs
>>>         ndvi.attrs['nodatavals'] = (0)
>>>         
>>>         # Write the data to a GeoTiff
>>>         ndvi.gw.to_raster('ndvi.tif',
>>>                           overwrite=True,
>>>                           n_workers=4,    # number of process workers sent to ``concurrent.futures``
>>>                           n_threads=2,    # number of thread workers sent to ``dask.compute``
>>>                           n_chunks=200,   # number of window chunks
>>>                           compress='lzw')

[jgrss]

The band math attributes below will now be updated.

new_attrs['nodatavals'] = (nodata)
new_attrs['scales'] = (1.0)
new_attrs['offsets'] = (0.0)

I also added a 'vi' attribute. When I merge your examples I will push these changes with version 1.2.5.

[jgrss]

I am closing this for now, but we can reopen if there are additional attributes that we want to add.