#!/bin/bash -x

#input is monthly GWP file
input=$1
indir=$2
outdir=$3

export GDAL_CACHEMAX=6000

date=$(basename $input .tif | cut -d'_' -f2)

yr=${date:0:4}
mn=${date:4:2}

yr_mn=${fn_date:0:6}

dir_to_search=${indir}/ocean_mask/
cp ${dir_to_search}/ocean_mask_global_GWP_res.tif .
ls .

if [ ! -f ${outdir}/GWP.OSWF.MONTHLY.${yr_mn}.v1.tif ];
then
	mkdir -p output/
	#reproject to WGS84
	gdalwarp -t_srs EPSG:4326 $input swat_${yr_mn}_global.tif -te -180 -90 180 90 -ts 172800 86400 -co compress=deflate -r near -srcnodata -3.4e+38 -dstnodata 255
	#change nodata value
	gdal_translate swat_${yr_mn}_global.tif swat_${yr_mn}_global_1.tif -a_nodata 254 -ot Byte -co compress=deflate
	
	#account for number of days/month
	if [[ $mn == "01" || $mn == "03" || $mn == "05" || $mn == "07" || $mn == "08" || $mn == "10" || $mn == "12" ]];
	then
		#change pixel value 255 to 0
		gdal_calc.py --outfile=swat_${yr_mn}_global_2.tif -A swat_${yr_mn}_global_1.tif --calc="0*(A==255)+A*logical_and(A>0, A<255)" --co=COMPRESS=LZW --type=Byte
		#add ocean mask to GWP file
		gdal_calc.py --outfile=swat_masked_${yr_mn}.tif -A swat_${yr_mn}_global_2.tif -B ocean_mask_global_GWP_res.tif --calc="A*B+31*(B==0)" --co=COMPRESS=LZW --type=Byte

	elif [ $mn == "02" ];
	then
		#account for leap years
		if [[ $yr_mn == "200002" || $yr_mn == "200402" || $yr_mn == "200802" || $yr_mn == "201202" || $yr_mn == "201602" || $yr_mn == "202002" || $yr_mn == "202402" ]];
		then
			#change pixel value 255 to 0
			gdal_calc.py --outfile=swat_${yr_mn}_global_2.tif -A swat_${yr_mn}_global_1.tif --calc="0*(A==255)+A*logical_and(A>0, A<255)" --co=COMPRESS=LZW --type=Byte
			#add ocean mask to GWP file
			gdal_calc.py --outfile=swat_masked_${yr_mn}.tif -A swat_${yr_mn}_global_2.tif -B ocean_mask_global_GWP_res.tif --calc="A*B+29*(B==0)" --co=COMPRESS=LZW --type=Byte
		else
			#change pixel value 255 to 0
			gdal_calc.py --outfile=swat_${yr_mn}_global_2.tif -A swat_${yr_mn}_global_1.tif --calc="0*(A==255)+A*logical_and(A>0, A<255)" --co=COMPRESS=LZW --type=Byte
			#add ocean mask to GWP file
			gdal_calc.py --outfile=swat_masked_${yr_mn}.tif -A swat_${yr_mn}_global_2.tif -B ocean_mask_global_GWP_res.tif --calc="A*B+28*(B==0)" --co=COMPRESS=LZW --type=Byte
		fi
	
	elif [[ $mn == "04" || $mn == "06" || $mn == "09" || $mn == "11" ]];
	then
		#change pixel value 255 to 0
		gdal_calc.py --outfile=swat_${yr_mn}_global_2.tif -A swat_${yr_mn}_global_1.tif --calc="0*(A==255)+A*logical_and(A>0, A<255)" --co=COMPRESS=LZW --type=Byte
		#add ocean mask to GWP file
		gdal_calc.py --outfile=swat_masked_${yr_mn}.tif -A swat_${yr_mn}_global_2.tif -B ocean_mask_global_GWP_res.tif --calc="A*B+30*(B==0)" --co=COMPRESS=LZW --type=Byte
	fi

	#change tif to cog file
	gdal_translate -of COG -co COMPRESS=ZSTD -co PREDICTOR=2 -co RESAMPLING=NEAREST swat_masked_${yr_mn}.tif GWP.OSWF.MONTHLY.${yr_mn}.v1.tif
	
	mv GWP.OSWF.MONTHLY.${yr_mn}.v1.tif output/
else
	echo "${outdir}/GWP.OSWF.MONTHLY.${yr_mn}.v1.tif already exists"
fi

exit 0