Skip to content

Commit

Permalink
Optimize planetry recipes for keeping highest resolution
Browse files Browse the repository at this point in the history
See #171 for background.  I have solved the issue by determining the exact increment and tile size required for the highest resolution, which gets a name based on the floor of the actual increment.  Need floor since 2.56 cannot be rint to 3 since we also will make a 03m grid.
Testing this now but seems to work so want to get these into the repo for further testing.
  • Loading branch information
PaulWessel committed Jan 9, 2023
1 parent 396da44 commit 96b90e0
Show file tree
Hide file tree
Showing 9 changed files with 66 additions and 25 deletions.
5 changes: 3 additions & 2 deletions recipes/mars_relief.recipe
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
# Recipe file for down-filtering MOLA
# 2023-01-07 PW
# 2023-01-09 PW
#
# We use a precision of 0.5 m with a 6000 m offset to fit the range of -8528 to +21226 in 16-bit ints
#
Expand Down Expand Up @@ -28,8 +28,9 @@
# DST_SRTM=no
#
# List of desired output resolution and chunk size. Flag the source resolution with code == master
# Given dimension 106694 x 53347 we have spacing of 12.1468873601 arc seconds as master and 14x7 tiles
# res unit tile chunk code
12 s 10 4096 master
12.1468873601 s 25.7142857143 4096 master
15 s 10 4096
30 s 15 4096
01 m 30 4096
Expand Down
3 changes: 2 additions & 1 deletion recipes/mercury_relief.recipe
Original file line number Diff line number Diff line change
Expand Up @@ -28,8 +28,9 @@
# DST_SRTM=no
#
# List of desired output resolution and chunk size. Flag the source resolution with code == master
# Given dimension 23040 x 11520 we have spacing of 56.25 arc seconds as master and 12x6 tiles
# res unit tile chunk code
56 s 30 4096 master
56.25 s 30 4096 master
01 m 30 4096
02 m 60 4096
03 m 90 2048
Expand Down
5 changes: 3 additions & 2 deletions recipes/moon_relief.recipe
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
# Recipe file for down-filtering LOLA
# 2023-01-03 PW
# 2023-01-09 PW
#
# We use a precision of 0.5 m with a zero offset to fit the range of -9128.5 to 10781.5 in 16-bit ints
#
Expand Down Expand Up @@ -28,8 +28,9 @@
# DST_SRTM=no
#
# List of desired output resolution and chunk size. Flag the source resolution with code == master
# Given dimension 92160 x 46080 we have spacing of 114.0625 arc seconds as master and 36x18 tiles
# res unit tile chunk code
14 s 10 4096 master
14.0625 s 10 4096 master
15 s 10 4096
30 s 15 4096
01 m 30 4096
Expand Down
3 changes: 2 additions & 1 deletion recipes/pluto_relief.recipe
Original file line number Diff line number Diff line change
Expand Up @@ -28,8 +28,9 @@
# DST_SRTM=no
#
# List of desired output resolution and chunk size. Flag the source resolution with code == master
# Given dimension 24888 x 12444 we have spacing of 52.07... arc seconds as master and 12x6 tiles
# res unit tile chunk code
52 s 15 4096 master
52.0732883317 s 30 4096 master
01 m 30 4096
02 m 60 4096
03 m 90 2048
Expand Down
12 changes: 10 additions & 2 deletions recipes/venus_relief.recipe
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
# Recipe file for down-filtering Magellan
# 2023-01-07 PW
# 2023-01-09 PW
#
# We use a precision of 0.5 m with a 4000 m offset to fit the range of -2951 to 11687 in 16-bit ints
#
Expand Down Expand Up @@ -28,7 +28,15 @@
# DST_SRTM=no
#
# List of desired output resolution and chunk size. Flag the source resolution with code == master
# Given dimension 8192 x 4096 we have spacing of 2.63... arc minutes as master and 6x4 tiles
# res unit tile chunk code
23 m 0 4096 master
2.63671875 m 60 4096 master
03 m 90 2048
04 m 180 2048
05 m 180 2048
06 m 0 4096
10 m 0 4096
15 m 0 4096
20 m 0 4096
30 m 0 4096
01 d 0 4096
17 changes: 10 additions & 7 deletions scripts/get_prefix_func.sh
Original file line number Diff line number Diff line change
@@ -1,17 +1,20 @@
#!/bin/bash
# Function that creates tile label
# Function that creates tile label in integer degrees W/S

function get_prefix () { # Takes west (in -180/180 range) and south and makes the {N|S}yy{W|E}xxx prefix
if [ $1 -ge 0 ]; then
X=$(printf "E%03d" $1)
# Get nearest integer degrees
W=$(gmt math -Q $1 RINT =)
S=$(gmt math -Q $2 RINT =)
if [ $W -ge 0 ]; then
X=$(printf "E%03d" $W)
else
t=$(gmt math -Q $1 NEG =)
t=$(gmt math -Q $W NEG =)
X=$(printf "W%03d" $t)
fi
if [ $2 -ge 0 ]; then
Y=$(printf "N%02d" $2)
if [ $S -ge 0 ]; then
Y=$(printf "N%02d" $S)
else
t=$(gmt math -Q $2 NEG =)
t=$(gmt math -Q $S NEG =)
Y=$(printf "S%02d" $t)
fi
echo ${Y}${X}
Expand Down
14 changes: 11 additions & 3 deletions scripts/srv_downsampler_grid.sh
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,13 @@
# title, radius of the planetary body, desired node registration and resolutions,
# desired output grid format and name prefix, and filter type, etc. Thus, this
# script should handle data from different planets.
# Note: If the highest resolution grid is not an integer unit then some exploration
# needs to be done to determine what increment and tile size give an integer number
# of tiles over 360 and 180 ranges. E.g., below is the master line for mars_relief
# (which had 200 m pixels on Mars spheroid) and earth_relief (which as 15s exactly):
# 12.1468873601 s 25.7142857143 4096 master
# 15 s 10 4096 master
# Easiest to work with number of rows and find suitable common factors.

if [ $# -eq 0 ]; then
echo "usage: srv_downsampler_grid.sh recipefile"
Expand Down Expand Up @@ -161,11 +168,12 @@ while read RES UNIT DST_TILE_SIZE CHUNK MASTER; do
echo "Bad resolution $RES - aborting"
exit -1
fi
if [ ! ${RES} = "01" ]; then # Use plural unit
IRES=$(gmt math -Q ${RES} FLOOR =)
if [ ${IRES} -gt 1 ]; then # Use plural unit
UNIT_NAME="${UNIT_NAME}s"
fi
for REG in ${DST_NODES}; do # Probably doing both pixel and gridline registered output, except for master */
DST_FILE=${DST_PLANET}/${DST_PREFIX}/${DST_PREFIX}_${RES}${UNIT}_${REG}.grd
DST_FILE=${DST_PLANET}/${DST_PREFIX}/${DST_PREFIX}_${IRES}${UNIT}_${REG}.grd
grdtitle="${TITLE} at ${RES} arc ${UNIT_NAME}"
# Note: The ${SRC_ORIG/+/\\+} below is to escape any plus-symbols in the file name with a backslash so grdedit -D will work
if [ -f ${DST_FILE} ]; then # Do nothing
Expand All @@ -177,7 +185,7 @@ while read RES UNIT DST_TILE_SIZE CHUNK MASTER; do
remark="Reformatted from master file ${SRC_ORIG/+/\\+} [${REMARK}]"
gmt grdedit ${DST_FILE} -D+t"${grdtitle}"+r"${remark}"+z"${SRC_NAME} (${SRC_UNIT})"
fi
elif [ "${UNIT}" = "s" ] && [ ${RES} -le 30 ]; then # Special handling of xxs -> 15s or 30s filter due to 64-bit bug in grdfilter?
elif [ "${UNIT}" = "s" ] && [ ${IRES} -le 30 ]; then # Special handling of xxs -> 15s or 30s filter due to 64-bit bug in grdfilter?
# See https://github.com/GenericMappingTools/remote-datasets/issues/32 - we do south and north hemisphere separately
# Get suitable Gaussian full-width filter rounded to nearest 0.01 km after adding 50 meters for noise
echo "Down-filter ${SRC_FILE} to ${DST_FILE}=${DST_MODIFY}"
Expand Down
12 changes: 10 additions & 2 deletions scripts/srv_downsampler_image.sh
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,13 @@
# format, radius of the planetary body, desired node registration and resolutions,
# name prefix, and filter type, etc. Thus, this script should handle images from
# different planets.
# Note: If the highest resolution image is not an integer unit then some exploration
# needs to be done to determine what increment and tile size give an integer number
# of tiles over 360 and 180 ranges. E.g., below is the master line for mars_relief
# (which had 200 m pixels on Mars spheroid) and earth_relief (which as 15s exactly):
# 12.1468873601 s 25.7142857143 4096 master
# 15 s 10 4096 master
# Easiest to work with number of rows and find suitable common factors.

if [ $# -eq 0 ]; then
echo "usage: srv_downsampler_image.sh recipefile"
Expand Down Expand Up @@ -102,10 +109,11 @@ while read RES UNIT TILE MASTER; do
echo "Bad resolution $RES - aborting"
exit -1
fi
if [ ! ${RES} = "01" ]; then # Use plural unit
IRES=$(gmt math -Q ${RES} FLOOR =)
if [ ${IRES} -gt 1 ]; then # Use plural unit
UNIT_NAME="${UNIT_NAME}s"
fi
DST_FILE=${DST_PLANET}/${DST_PREFIX}/${DST_PREFIX}_${RES}${UNIT}.tif
DST_FILE=${DST_PLANET}/${DST_PREFIX}/${DST_PREFIX}_${IRES}${UNIT}.tif
if [ -f ${DST_FILE} ]; then # Do nothing
echo "${DST_FILE} exist - skipping"
elif [ "X${MASTER}" = "Xmaster" ]; then # Just make a copy of the master to a new output file
Expand Down
20 changes: 15 additions & 5 deletions scripts/srv_tiler.sh
Original file line number Diff line number Diff line change
Expand Up @@ -16,6 +16,14 @@
# NOTE: We will ONLY look for the global files on this local machine. We first
# look in the staging/<planet> directory, and if not there then we look in the
# users server directory.
#
# Note: If the highest resolution grid is not an integer unit then some exploration
# needs to be done to determine what increment and tile size give an integer number
# of tiles over 360 and 180 ranges. E.g., below is the master line for mars_relief
# (which had 200 m pixels on Mars spheroid) and earth_relief (which as 15s exactly):
# 12.1468873601 s 25.7142857143 4096 master
# 15 s 10 4096 master
# Easiest to work with number of rows and find suitable common factors.

if [ $# -eq 0 ]; then
echo "usage: srv_tiler.sh recipe [-f]"
Expand Down Expand Up @@ -110,12 +118,13 @@ while read RES UNIT DST_TILE_SIZE CHUNK MASTER ; do
fi
for REG in ${DST_NODES}; do # Probably doing both pixel and gridline registered output, except for master */
# Name and path of grid we wish to tile
DST_TILE_TAG=${DST_PREFIX}_${RES}${UNIT}_${REG}
IRES=$(gmt math -Q ${RES} FLOOR =)
DST_TILE_TAG=${DST_PREFIX}_${IRES}${UNIT}_${REG}
DST_FILE=${DST_TILE_TAG}.grd
if [ -f ${DATADIR}/${DST_FILE} ]; then # found locally
DATAGRID=${DATADIR}/${DST_FILE}
else # Get it via local server files
DATAGRID=~/.gmt/server/${DST_PREFIX}_${RES}${UNIT}_${REG}.grd
DATAGRID=~/.gmt/server/${DST_PREFIX}_${IRES}${UNIT}_${REG}.grd
fi
if [ ! -f ${DATAGRID} ]; then # No
echo "No such file to tile: ${DATAGRID}"
Expand All @@ -126,7 +135,8 @@ while read RES UNIT DST_TILE_SIZE CHUNK MASTER ; do
FILTER_WIDTH=$(gmt math -Q ${SRC_RADIUS} 2 MUL PI MUL 360 DIV $INC MUL 0.05 ADD 10 MUL RINT 10 DIV =)
# Compute number of tiles required for this grid given nominal tile size.
# We enforce square tiles by only solving for ny and doubling it for nx
if [ $DST_TILE_SIZE -gt 0 ]; then # OK, we need to split the file into separate tiles
IDST_TILE_SIZE=$(gmt math -Q ${DST_TILE_SIZE} RINT =)
if [ ${IDST_TILE_SIZE} -gt 0 ]; then # OK, we need to split the file into separate tiles
ny=$(gmt math -Q 180 ${DST_TILE_SIZE} DIV =)
nx=$(gmt math -Q ${ny} 2 MUL =)
n_tiles=$(gmt math -Q $nx $ny MUL =)
Expand All @@ -135,7 +145,7 @@ while read RES UNIT DST_TILE_SIZE CHUNK MASTER ; do
# Build the list of w/e/s/n for the tiles
gmt grdinfo ${DATAGRID} -I${DST_TILE_SIZE} -D -C > ${TMP}/wesn.txt
# Determine local temporary tile directory for this product
TILEDIR=./${DST_PLANET}/${DST_PREFIX}/${DST_PREFIX}_${RES}${UNIT}_${REG}
TILEDIR=./${DST_PLANET}/${DST_PREFIX}/${DST_PREFIX}_${IRES}${UNIT}_${REG}
rm -rf ${TILEDIR}
mkdir -p ${TILEDIR}
while read w e s n; do
Expand All @@ -156,7 +166,7 @@ while read RES UNIT DST_TILE_SIZE CHUNK MASTER ; do
MSG="${TITLE} at ${RES}x${RES} arc ${UNAME} reduced by Gaussian ${DST_MODE} filtering (${FILTER_WIDTH} km fullwidth)"
fi
printf "/server/%s/%s/\t%s_%s_%s/\t%s\t%s\t%s\t%s\t%4s\t%s\t%s\t-\t-\t%s\t%s [%s]\n" \
${DST_PLANET} ${DST_PREFIX} ${DST_PREFIX} ${TAG} ${REG} ${TAG} ${REG} ${DST_SCALE} ${DST_OFFSET} ${SIZE} ${DST_TILE_SIZE} ${creation_date} ${DST_CPT} "${MSG}" "${CITE}" >> ${DST_PREFIX}_server.txt
${DST_PLANET} ${DST_PREFIX} ${DST_PREFIX} ${TAG} ${IREG} ${TAG} ${REG} ${DST_SCALE} ${DST_OFFSET} ${SIZE} ${DST_TILE_SIZE} ${creation_date} ${DST_CPT} "${MSG}" "${CITE}" >> ${DST_PREFIX}_server.txt

# Move the tiled grid away from this tree
mkdir -p ${TOPDIR}/staging/tiled
Expand Down

0 comments on commit 96b90e0

Please sign in to comment.