From ee2e30f6b35691394ca9b4ef8969a974578ae0d8 Mon Sep 17 00:00:00 2001 From: Quinn Hart Date: Tue, 1 Oct 2019 00:47:13 -0700 Subject: [PATCH] Updated goes for dynamic solar processing. Fixes #90 --- g.cimis/g.cimis/etc/configure.mk | 4 +- g.cimis/g.cimis/etc/goes.mk | 181 +++++++++++++++++++------ g.cimis/g.cimis/etc/solar_functions.mk | 1 + 3 files changed, 145 insertions(+), 41 deletions(-) diff --git a/g.cimis/g.cimis/etc/configure.mk b/g.cimis/g.cimis/etc/configure.mk index 453e0d2..83208a7 100755 --- a/g.cimis/g.cimis/etc/configure.mk +++ b/g.cimis/g.cimis/etc/configure.mk @@ -94,8 +94,8 @@ etc:=$(loc)/$(MAPSET)/etc SQLITE:=sqlite3 db.connect.database:=${loc}/${MAPSET}/sqlite.db -.PHONY: info -info:: +.PHONY: infox +infox:: @echo "#### g.gisenv: ####" @g.gisenv; @echo "#### g.gisenv: ####" diff --git a/g.cimis/g.cimis/etc/goes.mk b/g.cimis/g.cimis/etc/goes.mk index 883a6d3..fc01c5c 100755 --- a/g.cimis/g.cimis/etc/goes.mk +++ b/g.cimis/g.cimis/etc/goes.mk @@ -1,5 +1,88 @@ #! /usr/bin/make -f +define pod + +=pod + +=head1 SYNOPSIS + + g.cimis sec=goes cmd= [args=files=I] [arg=mapsets=I] + where command is one of: import, goes_to_solar, solar + +This g.cimis section (a Makefile) is used import GOES-1[67] data into a grass database, and potentially +project it into the solar mapset, and run solar calculations. It can work on either mapsets, or a list of files, or a +combination of both. It should work for multiple mapsets of files spanning multiple days. It will also work +specifing files one at a time, as in a incrontab. + +=head1 COMMANDS + +=over 4 + +=item B + +Show this information page + +=item B + +Import all specified I into the goes location. The files are imported into mapsets corresponding to each day, and +named at the time of acqusition. The I has no effect in this case. Also, files will be imported regardless of whether +they will be used in any subsequent solar calculation. That means nightime files will be included. + +=item B + +Project all rasters in the goes location into an equivalent structure in the I location. Any specified I will +first be imported, if they have not yet been. Both I and I parameters can be specified. + +=item B + +Compute all solar parameters for all rasters specified by either the I or I parameters. Files will first be +imported and projected into the solar space. The solar calculation expects that it knows all the files up to the latest one. +Running this command where some files/rasters are missing will cause errors in the calculation, as the algorythm needs a proper +order for the files. The I command can be modified with the B and B parameters as described +below. Note, that unlike the I command, the I command will not import night-time (unused) files. + +=back + +=head1 OPTIONS + +=over 4 + +=item <-n> + +This option is specified in the C. It will cause the operations that would be performed, to be written +to stdout. This is a good debugging tool. As each section in g.cimis is simply a makefile, this is equivalant +to C. + +=item B + +This is a space separated list of the GOES.pgm files that are ready for import into the system. Please note that this +program only imports pgm files that are in the very specific cookie-cutter regions used in cimis processing. These are +specified in the C file, in the extension. + +=item B + +This is a space separated list of mapsets within the B location. Every matching raster (glob pattern=????PST-B2) +will be included in the subsequent command. + +=item B + +Only used for the B command, if the I option is included, then all mapsets will be checked that a complete +days calculations have been performed, including the final sunset (ssetr-G) total radiance example. This is done, by checking +the latest raster time for each mapset, and adding a sunset calculation if they are all in the day. This should be added when +you want to calculate previously collected data. It should b be used in an incrontab type process. + +=item B + +B This parameter will only check for files in the B location, and not back to the B +location. This is useful if the original B data is somehow missing or removed. You can, and probably should use the +B parameter as well. + +=back + +=cut + +endef + ifndef configure.mk include configure.mk endif @@ -10,33 +93,38 @@ endif goes.mk:=1 + # Handy oneline functions f_fn=$(notdir $(basename $1)) f_mapset=$(word 1,$(subst T, ,$(call f_fn,$1))) f_rastname=$(word 2,$(subst T, ,$(call f_fn,$1)))T$(word 3,$(subst T, ,$(call f_fn,$1))) +f_rast_mapset=$(call f_rastname,$1)@$(call f_mapset,$1) -f_goes_rast=${GISDBASE}/${goes.loc}/$(call f_mapset,$1)/cellhd/$(call f_rastname,$1) -f_solar_rast=${GISDBASE}/${solar.loc}/$(call f_mapset,$1)/cellhd/$(call f_rastname,$1) #files:=$(wildcard /home/cimis/CA/*.pgm) -mapsets:=$(sort $(foreach f,${files},$(call f_mapset,$f))) +# you can call this two parameters; files, and/or mapsets. They will be combined together to create all rules +# Needed to import, project, and calculate solar parameters +file_mapsets:=$(sort $(foreach f,${files},$(call f_mapset,$f))) +all_mapsets:=$(sort ${mapsets} ${file_mapsets}) + + +info:: + @pod2usage -exit 0 $(firstword ${MAKEFILE_LIST}) + +#check:: +# @podchecker ${MAKEFILE_LIST} info:: - @echo goes.mk - @echo "Import GOES data into Grass" - @echo filenames:${filenames} - @echo mapset:${mapsets} - @echo $(call f_mapset,$(firstword ${files})) - @echo $(call f_rast,$(firstword ${files})) + @echo files:${files} + @echo mapsets:${all_mapsets} .PHONY: import solar goes_to_solar define _import -$(eval rast:=${GISDBASE}/${solar.loc}/$(call f_mapset,$1)/cellhd) -import::$(call f_goes_rast,$1) +import::${GISDBASE}/${goes.loc}/$(call f_mapset,$1)/cellhd/$(call f_rastname,$1) @echo "imported" -$(call f_goes_rast,$1):$1 +${GISDBASE}/${goes.loc}/$(call f_mapset,$1)/cellhd/$(call f_rastname,$1):$1 @$(call g.mapset-c,${goes.loc},$(call f_mapset,$1));\ echo -e '${import.wld}' > $(patsubst %.pgm,%.wld,$1);\ r.in.gdal --quiet --overwrite -o input=$1 output=$(call f_rastname,$1);\ @@ -46,46 +134,61 @@ $(call f_goes_rast,$1):$1 endef define _goes_to_solar -$(eval rast:=${GISDBASE}/${solar.loc}/$(call f_mapset,$1)/cellhd) +goes_to_solar::${GISDBASE}/${solar.loc}/$1/cellhd/$2 -goes_to_solar::$(call f_solar_rast,$1) - -$(call f_solar_rast,$1):$(call f_goes_rast,$1) - @$(call g.mapset-c,${solar.loc},$(call f_mapset,$1));\ +${GISDBASE}/${solar.loc}/$1/cellhd/$2:${GISDBASE}/${goes.loc}/$1/cellhd/$2 + @$(call g.mapset-c,${solar.loc},$1);\ g.region --quiet ${solar.region};\ - r.proj --quiet location=${goes.loc} mapset=$(call f_mapset,$1) \ - input=$(call f_rastname,$1) method=${solar.proj.method};\ - echo ${solar.loc}/$(call f_mapset,$1)/$(call f_rastname,$1) - + r.proj --quiet location=${goes.loc} mapset=$1 \ + input=$2 method=${solar.proj.method};\ + echo ${solar.loc}/$1/$2 endef -define _add_night -$(eval $(call _import,$1)) +define _add_predawn endef -define _add_day - -$(eval $(call _import,$1)) -$(eval $(call _goes_to_solar,$1)) -# This adds in 'solar' rule. Can't use with daily_solar as well, -# It would get overwritten. -$(eval $(call next_solar_calc,$(call f_rastname,$1),$(call f_mapset,$1))) +define _add_sunrise +$(eval $(call _goes_to_solar,$1,$2) $(call _sunrise,$2,$3,$1)) endef -define _add_sunrise -$(call _add_day,$1) +# This adds in 'solar' rule. Can't use with daily_solar as well, +# It would get overwritten. +# solar_functions go rast,prev,$mapset? +define _add_day +$(eval $(call _goes_to_solar,$1,$2) $(call _day,$2,$3,$1)) endef define _add_sunset -$(call _add_day,$1) +$(eval $(call _sunset,$2,$3,$1)) endef -# We have to be in the solar location to get proper sunrise time. -define import -$(eval tod:=$(shell $(call g.mapset-c,${solar.loc},$(call f_mapset,$1)); g.solar_time cmd=risedayset mapset=$(call f_mapset,$1) rast=$(call f_rastname,$1))) -$(call _add_$(firstword ${tod}),$1) +define _add_night endef -$(foreach m,$(sort $(foreach f,${files},$(call f_mapset,$f))),$(eval $(call mapset_targets,$m))) -$(foreach f,$(filter %.pgm,${files}),$(eval $(call import,$f))) +# Let's get all files that are either in the mapsets, or in the files, ready to be added. They are saved in a parameter list +# one for each mapset +$(foreach m,${all_mapsets},$(eval $m.files:=$(sort $(shell r.proj -l location=${goes.loc} mapset=$m 2>/dev/null || true | grep ????PST-B2) $(patsubst %@$m,%,$(filter %@$m,$(foreach i,$(filter %PST-B2.pgm,${files}),$(call f_rast_mapset,$i))))))) + +# And calculate the sunrise and sunset +$(foreach m,${all_mapsets},$(eval $m.sretr=$(shell m=$m; y=$${m%????}; md=$${m#????};m=$${md%??};d=$${md#??}; r.solpos -r year=$$y month=$$m day=$$d timezone=-8 | grep sretr_hhmm | cut -d= -f 2 | sed -e 's/://'))) +$(foreach m,${all_mapsets},$(eval $m.ssetr=$(shell m=$m; y=$${m%????}; md=$${m#????};m=$${md%??};d=$${md#??}; r.solpos -r year=$$y month=$$m day=$$d timezone=-8 | grep ssetr_hhmm | cut -d= -f 2 | sed -e 's/://'))) +#$(foreach m,${all_mapsets},$(info $m: sretr:${$m.sretr} ssetr:${$m.ssetr})) # Show them + +# And Calculate night/day for each file we are looking at +$(foreach m,${all_mapsets},$(eval $m.daynight:=$(foreach f,$(patsubst %PST-B2,%,${$m.files}),$(shell if [[ $f < ${$m.sretr} || $f = ${$m.sretr} ]]; then echo predawn; elif [[ "$f" < "${$m.ssetr}" ]]; then echo day; else echo night; fi)))) + +# And now find sunrise and sunset. Compare to the previous time / and it's type +$(foreach m,${all_mapsets},$(eval $m.prev:=0000 ${$m.files})) +$(foreach m,${all_mapsets},$(eval $m.prev_daynight:=x ${$m.daynight})) + +$(foreach m,${all_mapsets},$(eval $m.ranges:=$(foreach i,$(shell seq 1 $(words ${$m.files})),$(if $(filter predawn,$(word $i,${$m.daynight})),predawn:,$(if $(filter day,$(word $i,${$m.daynight})),$(if $(filter predawn,$(word $i,${$m.prev_daynight})),sunrise,day),$(if $(filter day,$(word $i,${$m.prev_daynight})),sunset,night))):$(word $i,${$m.files}):$(word $i,${$m.prev})))) + +# And finally create the proper rules, first mapset targets +$(foreach m,${all_mapsets},$(eval $(call mapset_targets,$m))) +# Then file imports +$(foreach f,$(sort $(filter %PST-B2.pgm,${files})),$(eval $(call _import,$f))) +# And also all the projection / and solar calculations +$(foreach m,${all_mapsets},$(foreach f,${$m.ranges},$(call _add_$(word 1,$(subst :, ,$f)),$m,$(word 2,$(subst :, ,$f)),$(word 3,$(subst :, ,$f))))) +# The user can provide a add-sunset parameter to put in a pseudo end of day if needed +$(if ${add-sunset},$(foreach m,${all_mapsets},$(if,$(filter-out night,$(lastword ${$m.daynight})),$(call _add_sunset,$m,${$m.ssetr},$(lastword ${$m.files}))))) diff --git a/g.cimis/g.cimis/etc/solar_functions.mk b/g.cimis/g.cimis/etc/solar_functions.mk index 1f8c0a7..bcbcc78 100644 --- a/g.cimis/g.cimis/etc/solar_functions.mk +++ b/g.cimis/g.cimis/etc/solar_functions.mk @@ -20,6 +20,7 @@ tl_DD:=01 01 01 01 07 07 07 07 07 07 \ 07 15 15 15 15 15 15 15 21 21 \ 21 21 21 21 28 28 28 28 28 28 28 + ###################################################### # # Daily Targets