[GIS] QGIS Temporal Controller with raster data

qgisqgis-3rasterspatio-temporal-datatemporal-controller

Starting from QGIS 3.14 TimeManager is dead, long live the Temporal Controller!
There are a lot of tutorials out there how "Temporal Controller" deals with vector layers, which looks nice. However, I'm working a lot with raster data and I cannot find any information, how to get a similar behaviour as it was in TimeManager.

I have a NetCDF like this:

ncdump -h SLA_ts.nc 
netcdf SLA_ts {
dimensions:
    time = 124 ;
    lat = 5 ;
    lon = 5 ;
variables:
    char crs ;
        crs:grid_mapping_name = "latitude_longitude" ;
        crs:long_name = "CRS definition" ;
        crs:semi_major_axis = 6378136.3 ;
        crs:inverse_flattening = 298.257 ;
    float time(time) ;
        time:long_name = "time" ;
        time:calendar = "standard" ;
        time:units = "days since 1858-11-17 00:00:00.0" ;
        time:actual_range = 52244., 55987.735625 ;
    float lat(lat) ;
        lat:standard_name = "latitude" ;
        lat:long_name = "latitude" ;
        lat:units = "degrees_north" ;
        lat:actual_range = 53.8, 54.2 ;
    float lon(lon) ;
        lon:standard_name = "longitude" ;
        lon:long_name = "longitude" ;
        lon:units = "degrees_east" ;
        lon:actual_range = 7.5, 7.9 ;
    float sla(time, lat, lon) ;
        sla:_FillValue = NaNf ;
        sla:grid_mapping = "crs" ;
        sla:long_name = "sea level anomaly" ;
        sla:unit = "m" ;

// global attributes:
        :Conventions = "CF-1.5" ;
}

When I load this file in QGIS, the time steps are loaded as bands (in the Symbology-Tab, they are displayed as Band 1: time=52244 (days since 1858-11-17 00:00:00) and so on…). However, the Temporal-Tab offers only a fixed time range as shown in the screenshot:enter image description here

I would expect that the Temporal Controller automatically recognizes that the layers are time steps but the button "Automatic" is inactive. Has anyone any experience with raster layers in "Temporal Controller"? How can I modify my data to be interpreted correctly? Or is this functionality simply not implemented (as the Temporal Controller is quite new)?

Best Answer

I wrote a PyQGIS script to handle the workaround for NetCDF files which uses temporal controller's updateTemporalRange signal. Just copy the following script to QGIS Python editor window, edit your dates, run the script and use Temporal Controller as usual.

import datetime

# These functions are part of https://github.com/GispoCoding/qgis_plugin_tools/blob/master/tools/raster_layers.py

def set_raster_renderer_to_singleband(layer: QgsRasterLayer, band: int = 1) -> None:
    """
    Set raster renderer to singleband
    :param layer: raster layer
    """
    # https://gis.stackexchange.com/a/377631/123927 and https://gis.stackexchange.com/a/157573/123927
    provider: QgsRasterDataProvider = layer.dataProvider()
    renderer: QgsSingleBandGrayRenderer = QgsSingleBandGrayRenderer(layer.dataProvider(), band)

    stats: QgsRasterBandStats = provider.bandStatistics(band, QgsRasterBandStats.All, layer.extent(), 0)
    min_val = max(stats.minimumValue, 0)
    max_val = max(stats.maximumValue, 0)

    enhancement = QgsContrastEnhancement(renderer.dataType(band))
    contrast_enhancement = QgsContrastEnhancement.StretchToMinimumMaximum
    enhancement.setContrastEnhancementAlgorithm(contrast_enhancement, True)
    enhancement.setMinimumValue(min_val)
    enhancement.setMaximumValue(max_val)
    layer.setRenderer(renderer)
    layer.renderer().setContrastEnhancement(enhancement)
    layer.triggerRepaint()
    
def set_band_based_on_range(layer: QgsRasterLayer, t_range: QgsDateTimeRange) -> int:
    """

    :param layer: raster layer
    :param t_range: temporal range
    :return: band number
    """
    band_num = 1
    tprops: QgsRasterLayerTemporalProperties = layer.temporalProperties()
    if tprops.isVisibleInTemporalRange(t_range) and t_range.begin().isValid() and t_range.end().isValid():
        if tprops.mode() == QgsRasterLayerTemporalProperties.ModeFixedTemporalRange:
            layer_t_range: QgsDateTimeRange = tprops.fixedTemporalRange()
            start: datetime.datetime = layer_t_range.begin().toPyDateTime()
            end: datetime.datetime = layer_t_range.end().toPyDateTime()
            delta = (end - start) / layer.bandCount()
            band_num = int((t_range.begin().toPyDateTime() - start) / delta) + 1
            set_raster_renderer_to_singleband(layer, band_num)
    return band_num
    
def set_fixed_temporal_range(layer: QgsRasterLayer, t_range: QgsDateTimeRange) -> None:
    """
    Set fixed temporal range for raster layer
    :param layer: raster layer
    :param t_range: fixed temporal range
    """
    mode = QgsRasterLayerTemporalProperties.ModeFixedTemporalRange
    tprops: QgsRasterLayerTemporalProperties = layer.temporalProperties()
    tprops.setMode(mode)
    if t_range.begin().timeSpec() == 0 or t_range.end().timeSpec() == 0:
        begin = t_range.begin()
        end = t_range.end()
        begin.setTimeSpec(Qt.TimeSpec(1))
        end.setTimeSpec(Qt.TimeSpec(1))
        t_range = QgsDateTimeRange(begin, end)
    tprops.setFixedTemporalRange(t_range)
    tprops.setIsActive(True)
    
def temporal_range_changed(t_range: QgsDateTimeRange):
    layer = iface.activeLayer()
    if isinstance(layer, QgsRasterLayer):
        set_band_based_on_range(layer, t_range)
    
def set_range():
    mode = QgsRasterLayerTemporalProperties.ModeFixedTemporalRange
    
temporal_controller: QgsTemporalController = iface.mapCanvas().temporalController()
temporal_controller.updateTemporalRange.connect(temporal_range_changed)
# Add one second to make the last frame visible
set_fixed_temporal_range(iface.activeLayer(), QgsDateTimeRange(datetime.datetime(2020, 11, 2, 15, 0), datetime.datetime(2020, 11, 3, 11, 0, 1))) 

Here is a simple gif animation showing the usage:

gif animation of usage

Related Question