Temporal data processing/seasonal aggregation

Seasonal aggregation example

This is an example showing the use of t.rast.aggregate.ds to aggregate a daily time series using a source time series with seasonal granularity. We create the daily maps for more than one year (400 maps). We then create the corresponding spatio-temporal raster data set (strds) and register maps in it.

# set the computational region and resolution
g.region -p n=-30 s=-50 e=-50 w=-70 res=1

# we do more than one year
for map in `seq 1 400` ; do
 # generate synthetic maps as a seeding
 r.mapcalc expression="daily_prec_${map} = 1"
 echo daily_prec_${map} >> map_list.txt
done

t.create type=strds temporaltype=absolute \
         output=precipitation_daily \
         title="Daily precipitation" \
         description="Test dataset with daily precipitation"

t.register -i type=raster input=precipitation_daily \
           file=map_list.txt start="2014-03-07" increment="1 days"

# check general information of the daily strds
t.info type=strds input=precipitation_daily

+-------------------- Space Time Raster Dataset -----------------------------+
|                                                                            |
+-------------------- Basic information -------------------------------------+
| Id: ........................ precipitation_daily@meteo
| Name: ...................... precipitation_daily
| Mapset: .................... meteo
| Creator: ................... veroandreo
| Temporal type: ............. absolute
| Creation time: ............. 2015-08-23 21:08:00.744761
| Modification time:.......... 2015-08-23 21:08:10.790248
| Semantic type:.............. mean
+-------------------- Absolute time -----------------------------------------+
| Start time:................. 2014-03-07 00:00:00
| End time:................... 2015-04-11 00:00:00
| Granularity:................ 1 day
| Temporal type of maps:...... interval
+-------------------- Spatial extent ----------------------------------------+
| North:...................... -30.0
| South:...................... -50.0
| East:.. .................... -50.0
| West:....................... -70.0
| Top:........................ 0.0
| Bottom:..................... 0.0
+-------------------- Metadata information ----------------------------------+
| Raster register table:...... raster_map_register_d14f30f1438840ac81cbfc0a554c1e69
| North-South resolution min:. 1.0
| North-South resolution max:. 1.0
| East-west resolution min:... 1.0
| East-west resolution max:... 1.0
| Minimum value min:.......... 1.0
| Minimum value max:.......... 1.0
| Maximum value min:.......... 1.0
| Maximum value max:.......... 1.0
| Aggregation type:........... None
| Number of registered maps:.. 400
|
| Title:
| Daily precipitation
| Description:
| Test dataset with daily precipitation
| Command history:
| # 2015-08-23 21:08:00 
| t.create type="strds" temporaltype="absolute"
|     output="precipitation_daily" title="Daily precipitation"
|     description="Test dataset with daily precipitation"
| # 2015-08-23 21:08:10 
| t.register -i type="raster"
|     input="precipitation_daily" file="map_list.txt" start="2014-03-07"
|     increment="1 days"
| 
+----------------------------------------------------------------------------+

Now, we'll create a second time series (might be raster, vector or raster3d) of seasonal granularity. We'll then take this granularity to aggregate our daily strds. Just as an example, we'll create a spatio-temporal vector data set (stvds).

# we create 4 maps of points (one per season)
for i in `seq 1 4` ; do 
 v.random output=points_${i} n=20
done

# create the stvds
t.create type=stvds temporaltype=absolute \
         output=points \
         title="Points" \
         description="Points for aggregation"

Now, we set the seasonal granularity for the vector maps. For this, we need to create a file with the respective map names and start and end date corresponding to astronomical seasons. We'll use this file as input for t.register.

points_1|2014-03-20|2014-06-21
points_2|2014-06-21|2014-09-23
points_3|2014-09-23|2014-12-21
points_4|2014-12-21|2015-03-20

t.register type=vector input=points file=seasonal_points.txt
t.info type=stvds input=points

+-------------------- Space Time Vector Dataset -----------------------------+
|                                                                            |
+-------------------- Basic information -------------------------------------+
| Id: ........................ points@meteo
| Name: ...................... points
| Mapset: .................... meteo
| Creator: ................... veroandreo
| Temporal type: ............. absolute
| Creation time: ............. 2015-08-23 21:10:11.462742
| Modification time:.......... 2015-08-23 21:49:59.554593
| Semantic type:.............. mean
+-------------------- Absolute time -----------------------------------------+
| Start time:................. 2014-03-20 00:00:00
| End time:................... 2015-03-20 00:00:00
| Granularity:................ 1 day
| Temporal type of maps:...... interval
+-------------------- Spatial extent ----------------------------------------+
| North:...................... -30.260367
| South:...................... -49.990657
| East:.. .................... -50.000035
| West:....................... -69.181185
| Top:........................ 0.0
| Bottom:..................... 0.0
+-------------------- Metadata information ----------------------------------+
| Vector register table:...... vector_map_register_727ad0a4157448ee8662d4345c248267
| Number of points ........... 80
| Number of lines ............ 0
| Number of boundaries ....... 0
| Number of centroids ........ 0
| Number of faces ............ 0
| Number of kernels .......... 0
| Number of primitives ....... 80
| Number of nodes ............ 0
| Number of areas ............ 0
| Number of islands .......... 0
| Number of holes ............ 0
| Number of volumes .......... 0
| Number of registered maps:.. 4
|
| Title:
| Points
| Description:
| Points for aggregation
| Command history:
| # 2015-08-23 21:10:11 
| t.create type="stvds" temporaltype="absolute"
|     output="points" title="Points" description="Points for aggregation"
| # 2015-08-23 21:49:59 
| t.register type="vector" input="points"
|     file="seasonal_points.txt"
| 
+----------------------------------------------------------------------------+

Finally, we aggregate our daily time series "precipitation daily" into seasonally accumulated precipitation with t.rast.aggregate.ds using method=sum.

# use everything between start date and end date
t.rast.aggregate.ds input=precipitation_daily \
                    output=precipitation_seasonal \
                    sample=points type=stvds \
                    base=prec_seasonal \
                    method=sum sampling=contains

# update information for the newly created seasonal time series
t.support input=precipitation_seasonal \
          title="Aggregated precipitation" \
          description="Aggregated seasonal precipitation dataset"

# check info                 
t.info type=strds input=precipitation_seasonal

+-------------------- Space Time Raster Dataset -----------------------------+
|                                                                            |
+-------------------- Basic information -------------------------------------+
| Id: ........................ precipitation_seasonal@meteo
| Name: ...................... precipitation_seasonal
| Mapset: .................... meteo
| Creator: ................... veroandreo
| Temporal type: ............. absolute
| Creation time: ............. 2015-08-23 21:50:26.134279
| Modification time:.......... 2015-08-23 21:50:34.365946
| Semantic type:.............. mean
+-------------------- Absolute time -----------------------------------------+
| Start time:................. 2014-03-20 00:00:00
| End time:................... 2015-03-20 00:00:00
| Granularity:................ 1 day
| Temporal type of maps:...... interval
+-------------------- Spatial extent ----------------------------------------+
| North:...................... -30.0
| South:...................... -50.0
| East:.. .................... -50.0
| West:....................... -70.0
| Top:........................ 0.0
| Bottom:..................... 0.0
+-------------------- Metadata information ----------------------------------+
| Raster register table:...... raster_map_register_be1ca044a59d49509e57a413e3781bc0
| North-South resolution min:. 1.0
| North-South resolution max:. 1.0
| East-west resolution min:... 1.0
| East-west resolution max:... 1.0
| Minimum value min:.......... 89.0
| Minimum value max:.......... 94.0
| Maximum value min:.......... 89.0
| Maximum value max:.......... 94.0
| Aggregation type:........... sum
| Number of registered maps:.. 4
|
| Title:
| Aggregated precipitation
| Description:
| Aggregated seasonal precipitation dataset
| Command history:
| # 2015-08-23 21:50:26 
| t.rast.aggregate.ds input="precipitation_daily"
|     output="precipitation_seasonal" sample="points" type="stvds"
|     base="prec_seasonal" method="sum" sampling="contains"
| # 2015-08-23 21:50:34 
| t.support input="precipitation_seasonal"
|     title="Aggregated precipitation"
|     description="Aggregated seasonal precipitation dataset"
| 
+----------------------------------------------------------------------------+

At the beginning we created 400 maps, remember??? We'll check now that only the maps "contained" in the period start_date-end_date of the stvds are used for the aggregation (the rest are left aside). That's why we use sampling=contains. Check t.rast.aggregate.ds for other method and sampling options.

# check that the sum is 365
t.rast.series in=precipitation_seasonal method=sum out=suma_prec_seasonal
d.mon wx0
d.rast suma_prec_seasonal
# use GUI tools to check raster values