Page MenuHomec4science
Files F111160102

compute_shading_stats.sh
No OneTemporary
Actions

File Metadata

Created
Tue, Apr 29, 23:01

compute_shading_stats.sh
View Options

#!/bin/bash
################## USER SETTINGS ##############
WORKDIR=$1
SOLPOS_FILE=$2
# Region_ID
REGION_ID=$3
# Layer names
ROOFTOP_MASK=$4
# Visibility parameters
MIN_VISIBILITY=$5
VIS_PERC=$6
## DIFFERENT FOR GVA VS CH!!
OUTPUT_BASE=$7
REGION_NAME=$8
#####
SUMMARY_STATS_BASE=$WORKDIR\/SHADE_$OUTPUT_BASE\_VIS_$VIS_PERC
MERGE_PYTHON_CODE=/home/walch/code/energy-potential/DOM/merge_shading_stats.py
echo Working in directory $WORKDIR
echo Reading from file $SOLPOS_FILE for region $REGION_ID
echo Mask layer: $ROOFTOP_MASK\, horizon basename: $OUTPUT_BASE\, region name: $REGION_NAME
echo Visibility parameters: $MIN_VISIBILITY\, $VIS_PERC
# ###################################################
## ------------ INPUT/OUTPUT FILES -------------##
# Solar positions:
LAYERFILE=$WORKDIR\/layers_$REGION_ID\.txt
echo Saving names of output layers in $LAYERFILE
# ############ Set UP GRASS ENVIRONMENT #############
mkdir $SUMMARY_STATS_BASE
## Set region
g.region region=$REGION_NAME
REGION=current
## Add mapset shading_rooftops
g.mapsets shade_rooftops operation=add
g.mapsets -p
#################################################
echo STARTING AT $(date)
# file that is used to store only the names of the layers, without month and hour tags
# this file is used as input to r.series
TMP_SOLPOS_FILE=$WORKDIR\/TMP_SOLPOS_$REGION_ID\.txt
cp $SOLPOS_FILE $TMP_SOLPOS_FILE
TMP_LAYER_FILE=$WORKDIR\/TMP_layers_$REGION_ID\.txt
# remove temporary files (in case they exist already)
rm $LAYERFILE
rm $TMP_LAYER_FILE
## for-loop from list: loop over all time stamps and create boolean masks indicating which pixels are visible
while read MONTH HOUR DIRECTION ALTITUDE DIRECTION_RND; do
echo Month: $MONTH\; Hour: $HOUR\; Direction: $DIRECTION\; Altitude: $ALTITUDE
## Set up variables and commands
HORIZON_RAW=RAW_$OUTPUT_BASE
HORIZON_RAW_OUT=$HORIZON_RAW\_$DIRECTION_RND\_0_ROOFS
HORIZON_MASKED_BOOL=$OUTPUT_BASE\_$MONTH\_$HOUR\h
EXPRESSION="$HORIZON_MASKED_BOOL = ($ALTITUDE>$HORIZON_RAW_OUT)"
## Execute
r.mapcalc expression="$EXPRESSION" region=$REGION
echo Computed visible points and applied mask. Saved in layer $HORIZON_MASKED_BOOL
echo $HORIZON_MASKED_BOOL >> $TMP_LAYER_FILE
echo $MONTH $HOUR $HORIZON_MASKED_BOOL >> $LAYERFILE
echo FINISHED ITERATION AT $(date)
done < $TMP_SOLPOS_FILE
echo FINISHED MASKING HORIZONS AT $(date)
echo
#########################################################
# Sum together all layers created in the previous loop
FILE_OF_SUMS=$OUTPUT_BASE\_MEAN
r.series -nz file=$TMP_LAYER_FILE output=$FILE_OF_SUMS method=average
rm $TMP_LAYER_FILE
# copy file to edit null values (set all values below visibility threshold to null)
TMP_FILE_OF_SUMS=TMP_$FILE_OF_SUMS
g.copy raster=$FILE_OF_SUMS\,$TMP_FILE_OF_SUMS
r.null $TMP_FILE_OF_SUMS setnull=0-$MIN_VISIBILITY
# compute summary statistics: null_cells are 'fully shaded', non_null_cells are partly shaded
r.univar -t $TMP_FILE_OF_SUMS zones=$ROOFTOP_MASK output=$SUMMARY_STATS_BASE\/summary_stats.csv --overwrite
# calculate a new mask of rooftops, with only not fully shaded cells
UNSHADED_MASK=$ROOFTOP_MASK\_vis_$VIS_PERC\_unshaded
EXPRESSION="$UNSHADED_MASK = int($ROOFTOP_MASK*($TMP_FILE_OF_SUMS*0+1))"
r.mapcalc expression="$EXPRESSION" region=$REGION --overwrite
echo FINISHED COMPUTING FULLY SHADED CELLS AT $(date)
echo
#########################################################
# compute statistics for all partly shaded roofs
## for-loop from list:
while read MONTH HOUR LAYER; do
echo Month: $MONTH\; Hour: $HOUR\; Layer: $LAYER\
r.univar -t map=$LAYER zones=$UNSHADED_MASK output=$SUMMARY_STATS_BASE\/stats_$MONTH\_$HOUR\h.csv --overwrite
done < $LAYERFILE
g.remove -f type=raster pattern=TMP_*
rm $LAYERFILE
rm $TMP_SOLPOS_FILE
echo FINISHED COMPUTING PARTLY SHADED CELLS AT $(date)
## RUN PYTHON SCRIPT merge_summary_files
source activate py2
python $MERGE_PYTHON_CODE $SUMMARY_STATS_BASE
source deactivate py2
echo FINISHED MERGING STATISTICS AT $(date)

Event Timeline

c4science · Help