|
22 | 22 | #
|
23 | 23 | ##########################################################################
|
24 | 24 |
|
| 25 | +# functions/definitions ------------------------------------------------------------------ |
| 26 | +export PROG=`basename $0`; |
| 27 | +export BIN="$( cd "$( dirname "${BASH_SOURCE[0]}" )" >/dev/null 2>&1 && pwd )" |
| 28 | + |
| 29 | +MANDATORY_ARGS=2 |
| 30 | + |
| 31 | +export CALC_EXE="gdal_calc.py" |
| 32 | +export PARALLEL_EXE="parallel" |
| 33 | + |
| 34 | +echoerr(){ echo "$PROG: $@" 1>&2; } # warnings and/or errormessages go to STDERR |
| 35 | +export -f echoerr |
| 36 | + |
| 37 | +export DEBUG=false # display debug messages? |
| 38 | +debug(){ if [ "$DEBUG" == "true" ]; then echo "DEBUG: $@"; fi } # debug message |
| 39 | +export -f debug |
| 40 | + |
| 41 | +cmd_not_found(){ # check required external commands |
| 42 | + for cmd in "$@"; do |
| 43 | + stat=`which $cmd` |
| 44 | + if [ $? != 0 ] ; then echoerr "\"$cmd\": external command not found, terminating..."; exit 1; fi |
| 45 | + done |
| 46 | +} |
| 47 | +export -f cmd_not_found |
| 48 | + |
| 49 | + |
| 50 | +help(){ |
| 51 | +cat <<HELP |
| 52 | +
|
| 53 | +Usage: $PROG [-sldobj] input-basename calc-expr |
| 54 | +
|
| 55 | + optional: |
| 56 | + -s = pixel resolution of cubed data, defaults to 10 |
| 57 | + -l = input-layer: band number in case of multi-band input rasters, |
| 58 | + defaults to 1 |
| 59 | + -d = input directory: the datacube directory |
| 60 | + defaults to current directory |
| 61 | + 'datacube-definition.prj' needs to exist in there |
| 62 | + -o = output directory: the directory where you want to store the cubes |
| 63 | + defaults to current directory |
| 64 | + -b = basename of output file (without extension) |
| 65 | + defaults to the basename of the input-file, |
| 66 | + appended by '_procmask' |
| 67 | + -j = number of jobs, defaults to 'as many as possible' |
| 68 | +
|
| 69 | + Positional arguments: |
| 70 | + - input-basename: basename of input data |
| 71 | + - calc-expr: Calculation in gdalnumeric syntax, e.g. 'A>2500'" |
| 72 | + The input variable is 'A' |
| 73 | + For details about GDAL expressions, see |
| 74 | + https://gdal.org/programs/gdal_calc.html |
| 75 | +
|
| 76 | + ----- |
| 77 | + see https://force-eo.readthedocs.io/en/latest/components/auxilliary/procmask.html |
| 78 | +
|
| 79 | +HELP |
| 80 | +exit 1 |
| 81 | +} |
| 82 | +export -f help |
| 83 | + |
| 84 | +# important, check required commands !!! dies on missing |
| 85 | +cmd_not_found "$CALC_EXE $PARALLEL_EXE" |
25 | 86 |
|
26 | 87 | function mask(){
|
27 | 88 |
|
28 | 89 | FILE=$1
|
29 | 90 | TILE=$(dirname $FILE)
|
30 | 91 |
|
31 |
| - mkdir -p $OUT/$TILE |
32 |
| - |
33 |
| - gdal_calc.py -A $FILE --A_band=$IBAND --outfile=$OUT/$TILE/$OBASE --calc=$CALC --NoDataValue=255 --type=Byte --format='GTiff' --creation-option='COMPRESS=LZW' --creation-option='PREDICTOR=2' --creation-option='NUM_THREADS=ALL_CPUS' --creation-option='BIGTIFF=YES' --creation-option="BLOCKXSIZE=$XBLOCK" --creation-option="BLOCKYSIZE=$YBLOCK" |
| 92 | + mkdir -p $DOUT/$TILE |
| 93 | + |
| 94 | + gdal_calc.py \ |
| 95 | + -A $FILE \ |
| 96 | + --A_band=$LAYER \ |
| 97 | + --outfile=$DOUT/$TILE/$FOUT \ |
| 98 | + --calc=$CALC \ |
| 99 | + --NoDataValue=255 \ |
| 100 | + --type=Byte \ |
| 101 | + --format='GTiff' \ |
| 102 | + --creation-option='COMPRESS=LZW' \ |
| 103 | + --creation-option='PREDICTOR=2' \ |
| 104 | + --creation-option='NUM_THREADS=ALL_CPUS' \ |
| 105 | + --creation-option='BIGTIFF=YES' \ |
| 106 | + --creation-option="BLOCKXSIZE=$XBLOCK" \ |
| 107 | + --creation-option="BLOCKYSIZE=$YBLOCK" |
34 | 108 |
|
35 | 109 | }
|
36 |
| - |
37 | 110 | export -f mask
|
38 | 111 |
|
39 | 112 |
|
40 |
| -EXPECTED_ARGS=7 |
41 |
| - |
42 |
| -if [ $# -ne $EXPECTED_ARGS ] |
43 |
| -then |
44 |
| - echo "Usage: `basename $0` input-cube output-cube" |
45 |
| - echo " input-base output-base" |
46 |
| - echo " input-band calc-expr resolution" |
47 |
| - echo "" |
48 |
| - echo " input-cube: directory of input cubes" |
49 |
| - echo " output-cube: directory of output cubes" |
50 |
| - echo " input-base: basename of cubed input raster" |
51 |
| - echo " output-base: basename of cubed processing masks" |
52 |
| - echo " input-base: band of cubed input raster, " |
53 |
| - echo " from which to generate the processing mask" |
54 |
| - echo " calc-expr: Calculation in gdalnumeric syntax, e.g. 'A>2500'" |
55 |
| - echo " The input variable is 'A'" |
56 |
| - echo " For details about GDAL expressions, see " |
57 |
| - echo " https://gdal.org/programs/gdal_calc.html" |
58 |
| - echo " resolution: the resolution of the cubed data" |
59 |
| - echo "" |
60 |
| - exit |
| 113 | +# now get the options -------------------------------------------------------------------- |
| 114 | +ARGS=`getopt -o hvis:l:d:o:b:j: --long help,version,info,resolution:,layer:,input:,output:,basename:,jobs: -n "$0" -- "$@"` |
| 115 | +if [ $? != 0 ] ; then help; fi |
| 116 | +eval set -- "$ARGS" |
| 117 | + |
| 118 | +RES=10 |
| 119 | +LAYER=1 |
| 120 | +DINP=$PWD |
| 121 | +DOUT=$PWD |
| 122 | +OBASE="DEFAULT" |
| 123 | +NJOB=0 |
| 124 | + |
| 125 | +while :; do |
| 126 | + case "$1" in |
| 127 | + -h|--help) help ;; |
| 128 | + -v|--version) echo "version-print to be implemented"; exit 0;; |
| 129 | + -i|--info) echo "Processing masks from raster images"; exit 0;; |
| 130 | + -s|--resolution) RES="$2"; shift ;; |
| 131 | + -l|--layer) LAYER="$2"; shift;; |
| 132 | + -d|--input) DINP=$(readlink -f" $2"); shift ;; |
| 133 | + -o|--output) DOUT=$(readlink -f "$2"); shift ;; |
| 134 | + -b|--basename) OBASE="$2"; shift ;; |
| 135 | + -j|--jobs) NJOB="$2"; shift ;; |
| 136 | + -- ) shift; break ;; |
| 137 | + * ) break ;; |
| 138 | + esac |
| 139 | + shift |
| 140 | +done |
| 141 | + |
| 142 | +if [ $# -ne $MANDATORY_ARGS ] ; then |
| 143 | + echoerr "Mandatory argument is missing."; help |
| 144 | +else |
| 145 | + IBASE="$1" # basename (with extension) |
| 146 | + CALC="$2" # GDAL calc expression |
61 | 147 | fi
|
62 | 148 |
|
63 |
| - |
64 |
| -INP=$1 |
65 |
| -OUT=$2 |
66 |
| -IBASE=$3 |
67 |
| -OBASE=$4 |
68 |
| -IBAND=$5 |
69 |
| -CALC=$6 |
70 |
| -RES=$7 |
71 |
| - |
72 |
| - |
73 |
| -if [ ! -r $INP ]; then |
74 |
| - echo "$INP is not existing/readable" |
75 |
| - exit |
| 149 | +debug "input dir: $DINP" |
| 150 | +debug "output dir: $DOUT" |
| 151 | +debug "input base: $IBASE" |
| 152 | +debug "output base: $OBASE" |
| 153 | +debug "input layer: $LAYER" |
| 154 | +debug "output res: $RES" |
| 155 | +debug "expression: $CALC" |
| 156 | +debug "njobs: $NJOB" |
| 157 | + |
| 158 | +if [ ! -r $DINP ]; then |
| 159 | + echoerr "$DINP is not existing/readable"; help |
76 | 160 | fi
|
77 | 161 |
|
78 |
| -if [ ! -w $OUT ]; then |
79 |
| - echo "$INP is not existing/writeable" |
80 |
| - exit |
| 162 | +if [ ! -w $DOUT ]; then |
| 163 | + echoerr "$DOUT is not existing/writeable"; help |
81 | 164 | fi
|
82 | 165 |
|
83 |
| -if [ ! -r $INP/datacube-definition.prj ]; then |
84 |
| - echo "$OUT/datacube-definition.prj is not existing/readable" |
85 |
| - exit |
| 166 | +if [ ! -r $DINP/datacube-definition.prj ]; then |
| 167 | + echoerr "$DINP/datacube-definition.prj is not existing/readable"; help |
86 | 168 | fi
|
87 | 169 |
|
88 |
| -if [ ! -r $OUT/datacube-definition.prj ]; then |
| 170 | +if [ -r $DOUT/datacube-definition.prj ]; then |
| 171 | + DIFF=$(diff "$DINP/datacube-definition.prj" "$DOUT/datacube-definition.prj") |
| 172 | + NUM=$(echo $DIFF | wc -w) |
| 173 | + if [ $NUM -gt 0 ]; then |
| 174 | + echoerr "input and output datacubes do not match"; help |
| 175 | + fi |
| 176 | +else |
89 | 177 | echo "copying datacube-definition.prj"
|
90 |
| - cp $INP/datacube-definition.prj $OUT/datacube-definition.prj |
| 178 | + cp "$DINP/datacube-definition.prj" "$DOUT/datacube-definition.prj" |
| 179 | +fi |
| 180 | + |
| 181 | +if [ $OBASE == "DEFAULT" ]; then |
| 182 | + CORE=${IBASE%%.*} # corename from input file |
| 183 | + FOUT="$CORE""_procmask.tif" |
| 184 | +else |
| 185 | + CORE=${OBASE%%.*} # corename from user parameters |
| 186 | + FOUT="$CORE.tif" |
91 | 187 | fi
|
92 | 188 |
|
93 |
| -# force tif extension |
94 |
| -OBASE=$(basename $OBASE) |
95 |
| -OBASE=${OBASE%%.*} |
96 |
| -OBASE=$OBASE".tif" |
| 189 | +debug "output file: $FOUT" |
97 | 190 |
|
| 191 | +# main thing ----------------------------------------------------------------------------- |
98 | 192 |
|
99 | 193 | NOW=$PWD
|
100 |
| -cd $INP |
| 194 | +cd $DINP |
101 | 195 |
|
102 | 196 | # list with input images
|
103 |
| -TEMP=temp-force-procmask.txt |
| 197 | +TEMP="$DOUT/temp-force-procmask.txt" |
104 | 198 |
|
105 |
| -ls X*/$IBASE > $TEMP |
| 199 | +ls X*/$IBASE > "$TEMP" |
106 | 200 |
|
107 |
| -if [ $(cat $TEMP | wc -l) -lt 1 ]; then |
108 |
| - echo "could not find any instance of $IBASE in $INP" |
109 |
| - rm $TEMP |
110 |
| - exit |
| 201 | +if [ $(cat "$TEMP" | wc -l) -lt 1 ]; then |
| 202 | + echoerr "could not find any instance of $IBASE in $DINP" |
| 203 | + rm "$TEMP" |
| 204 | + help |
111 | 205 | fi
|
112 | 206 |
|
113 | 207 | # tile /chunk size
|
114 |
| -TILESIZE=$(head -n 6 $OUT/datacube-definition.prj | tail -1 ) |
115 |
| -CHUNKSIZE=$(head -n 7 $OUT/datacube-definition.prj | tail -1 ) |
| 208 | +TILESIZE=$(head -n 6 $DOUT/datacube-definition.prj | tail -1 ) |
| 209 | +CHUNKSIZE=$(head -n 7 $DOUT/datacube-definition.prj | tail -1 ) |
116 | 210 |
|
117 | 211 | # block size
|
118 | 212 | XBLOCK=$(echo $TILESIZE $RES | awk '{print int($1/$2)}')
|
119 | 213 | YBLOCK=$(echo $CHUNKSIZE $RES | awk '{print int($1/$2)}')
|
120 | 214 |
|
121 |
| -export OUT=$OUT |
122 |
| -export OBASE=$OBASE |
123 |
| -export IBAND=$IBAND |
| 215 | +export DOUT=$DOUT |
| 216 | +export FOUT=$FOUT |
| 217 | +export LAYER=$LAYER |
124 | 218 | export CALC=$CALC
|
125 | 219 | export XBLOCK=$XBLOCK
|
126 | 220 | export YBLOCK=$YBLOCK
|
127 | 221 |
|
128 |
| -parallel -a $TEMP --eta mask {} |
| 222 | +$PARALLEL_EXE -j $NJOB -a $TEMP --eta mask {} |
129 | 223 |
|
130 | 224 | rm $TEMP
|
131 | 225 |
|
132 | 226 | cd $PWD
|
133 | 227 |
|
| 228 | +exit 0 |
0 commit comments