gmt-template/gmtxy-image2.sh
2020-04-26 12:09:20 +08:00

296 lines
9.8 KiB
Bash
Executable File
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

#!/bin/bash
# 包含GMT自带脚本文件 其中包含了一些有用的功能 比如获取网格文件的范围
. gmt_shell_functions.sh
# 包含dispOption脚本
. dispOptions
# 定义一个函数 执行输入的语句或者将其显示在屏幕上
RunOrEcho()
{
first_str=`echo ${2// /''}`
if [[ x${3} != x ]]; then
sec_str=`echo ${3// /''}`
fi
if [[ ${1} == 1 ]]; then
${first_str}
if [[ x${3} != x ]]; then
${sec_str}
fi
else
if [[ x${3} != x ]]; then
printf "%s\n%s\n" "${first_str}" "${sec_str}"
else
printf "%s\n" "${first_str}"
fi
fi
}
# GMT显示平面数据脚本输入文件为网格文件没有包含网格化语句因为网格化过程中的情况多样化建议在其他脚本中个别添加再调用此脚本
# 准备几个常用的排版参数组 记录的参数包括 图片的宽度 默认坐标标示设定 色标位置 色标的大小 标示字体大小
com_layouts=("1.5i,WesNZ,0.1i/-0.2i,1.3i/0.05i+h,10.5p" \
"1.5i,WeSnZ,1.65i/0.05i,1.3i/0.05i,10.5p" \
"1.2i,WesNZ,0.1i/-0.1i,1.0i/0.03i+h,7.5p" \
"1.2i,WeSnZ,1.3i/0.05i,1.0i/0.05i,7.5p")
# 初始化参数
data='null'
outpsfile='null'
gridData='null'
unit='m'
color='rainbow'
range='null'
overwriteRange='null'
overwrite=0
labels=("x (m)" "y (m)")
plotgrad=0
axistick=("a" "a") #横纵坐标轴间隔 a表示自动
bartick="a" #色标轴标轴间隔 a表示自动
frameset="null"
plot_colorbar=1
annotation="null"
keep_open=0
from_open=">"
image_origin=("0.5" "0.5")
bar_origin=("0" "0")
polyfile='null'
lay_num=0
run_command=1
# 从命令行获取参数
while getopts "hi:o:L:r:u:c:l:G:t:v:f:a:p:X:x:gnKOS" arg
do
case $arg in
h)
dispTitle "${0##*/}" "simple grid file mapping using GMT script under x-y plane. \
The template accepts a grid (.nc or .grd) file as input and outputs a .png and a .eps file. \
For further explanations, please look for GMT's manuscripts."
dispAuthorInfo "Dr. Yi Zhang (yizhang-geo@zju.edu.cn)"
dispUsage "${0##*/} -i<grid-data> [-o<img-file>] [-L<number>] [-r<xmin>/<xmax>/<ymin>/<ymax>] [-u<unit>] \
[-c<cpt-file>] [-t<number>,<number>] [-v<number>] [-l<x-label>;<y-label>] [-g] [-G<grad-file>] \
[-f<wEsNZ>] [-p<poly-file>] [-n] [-a<x-coor>,<y-coor>,<annotation>] [-K] [-O] [-S] \
[-X<x-offset>,<y-offset>] [-x<x-offset>,<y-offset>]"
dispOptionShort "-i" "Input grid file."
dispOptionShort "-o" "Output image file. File extension is not needed. the input filename will be \
used if this option is not set."
dispOptionShort "-L" "Select layouts. enter a number to select the layout to use (0-3). \
The default is zero."
dispOptionShort "-r" "Data range. The template will ditect the input data range automaticly, \
the use of this option will overwrite the range."
dispOptionShort "-u" "Data unit. The default is meter."
dispOptionShort "-c" "Color cpt. The default is rainbow."
dispOptionShort "-t" "Intervals of axis's ticks. the script will set the intervals automatically \
if -t option is not set."
dispOptionShort "-v" "Intervals of color bar's labels. the script will set the intervals \
automatically if -v option is not set."
dispOptionShort "-l" "Axis labels separated by semicolons. The deafults are x (m) and y (m)."
dispOptionShort "-g" "Use grdgradient. The default is false."
dispOptionShort "-G" "Provide a different input grid data for the use of grdgradient. \
This should be used with the '-g' option at the same time."
dispOptionShort "-f" "Set frames for mapping axises. The default is wEsNZ."
dispOptionShort "-n" "Don't plot color bar."
dispOptionShort "-a" "Add an annotation on the map. you need to enter a string contains both \
coordinates and text of the annotation, such as '10,10,(a)'."
dispOptionShort "-p" "Plot polygons via a file."
dispOptionShort "-K" "Keep the .ps file open."
dispOptionShort "-O" "Continue from a previous .ps file."
dispOptionShort "-X" "Move the starting point for mapping image."
dispOptionShort "-x" "Move the starting point for mapping color bar."
dispOptionShort "-S" "Instead of running the script, print a executable script on screen. \
This option is designed to help building data specified scripts using the template."
exit 0;;
i)
data=$OPTARG;;
o)
outpsfile=$OPTARG;;
L)
lay_num=$OPTARG;;
u)
unit=$OPTARG;;
c)
color=$OPTARG;;
t)
axistick=(${OPTARG//,/ });;
X)
image_origin=(${OPTARG//,/ });;
x)
bar_origin=(${OPTARG//,/ });;
v)
bartick=$OPTARG;;
f)
frameset=$OPTARG;;
r)
overwrite=1
overwriteRange=$OPTARG;;
l)
OLD_IFS="${IFS}"
IFS=";"
labels=(${OPTARG})
IFS="${OLD_IFS}";;
g)
plotgrad=1;;
n)
plot_colorbar=0;;
G)
gridData=$OPTARG;;
a)
annotation=$OPTARG;;
p)
polyfile=$OPTARG;;
K)
keep_open=1;;
O)
from_open="-O >>";;
S)
run_command=0;;
?)
printf "error: unknow argument\nuse -h option to see help information\n"
exit 1;;
esac
done
# 进行必要的参数检查
if [[ $data == "null" ]]; then
printf "error: no input file name\nuse -h option to see help information\n"
exit 1
else
# 初始化临时文件名
cptfile=user.cpt
if [[ ${outpsfile} == 'null' ]]; then
psfile=${data%.*}.ps
jpgfile=${data%.*}.png
else
psfile=${outpsfile}.ps
jpgfile=${outpsfile}.png
fi
# 获取排版
layouts=(${com_layouts[${lay_num}]//,/ })
# 获取坐标轴设置
if [[ ${frameset} == 'null' ]]; then
frameset=${layouts[1]}
fi
# 获取网格范围
if [[ $overwrite == 1 ]]; then
range=${overwriteRange}
else
range=$(gmt_get_gridregion ${data})
fi
#根据横纵坐标的范围计算图片的长和高 我们默认图片宽度为1.5i
range_coor=(${range//// })
pic_height=`echo "scale=4; ${layouts[0]%i} * (${range_coor[3]}-${range_coor[2]})/(${range_coor[1]}-${range_coor[0]})"|bc`
# 输出一个bash脚本的头部到屏幕
if [[ ${run_command} == 0 && ${from_open} == '>' ]]; then
echo "#!/bin/bash"
fi
# 设置绘图参数
if [[ ${from_open} == '>' ]]; then
RunOrEcho ${run_command} "gmt gmtset \
FONT_ANNOT_PRIMARY=${layouts[4]},Times-Roman,black \
MAP_FRAME_PEN=thinnest,black \
MAP_GRID_PEN_PRIMARY=thinnest,black \
MAP_TICK_PEN_PRIMARY=thinnest,black \
MAP_TICK_LENGTH_PRIMARY=1p/0.5p \
MAP_TITLE_OFFSET=7.5p \
MAP_GRID_CROSS_SIZE_PRIMARY=2p \
FONT_LABEL=${layouts[4]},Times-Roman,black \
MAP_FRAME_AXES=${frameset} \
MAP_LABEL_OFFSET=2.5p \
MAP_ANNOT_OFFSET_PRIMARY=2.5p"
fi
# 处理-B选项的参数
if [[ ${axistick[0]} == 'a' ]]; then
axistick[0]=${axistick[0]}g
axistick[1]=${axistick[1]}g
else
axistick[0]=${axistick[0]}g${axistick[0]}
axistick[1]=${axistick[1]}g${axistick[1]}
fi
RunOrEcho ${run_command} "gmt grd2cpt ${data} -C${color} -R${range} -Z -D > ${cptfile}"
if [[ $plotgrad == 1 ]]; then
gradfile=${data%.*}Grad.nc
if [[ $gridData == 'null' ]]; then
RunOrEcho ${run_command} "gmt grdgradient ${data} -G${gradfile} -Nt -A0/45"
else
RunOrEcho ${run_command} "gmt grdgradient ${gridData} -G${gradfile} -Nt -A0/45"
fi
if [[ ${run_command} == 0 ]]; then
echo "gmt grdimage ${data} -R${range} -C${cptfile} -I${gradfile} \
-Bx${axistick[0]}+l\"${labels[0]}\" -By${axistick[1]}+l\"${labels[1]}\" \
-JX${layouts[0]}/${pic_height}i -X${image_origin[0]}i -Y${image_origin[1]}i -K -P ${from_open} $psfile"
else
gmt grdimage ${data} -R${range} -C${cptfile} -I${gradfile} \
-Bx${axistick[0]}+l"${labels[0]}" -By${axistick[1]}+l"${labels[1]}" \
-JX${layouts[0]}/${pic_height}i -X${image_origin[0]}i -Y${image_origin[1]}i -K -P ${from_open} $psfile
fi
else
if [[ ${run_command} == 0 ]]; then
echo "gmt grdimage ${data} -R${range} -C${cptfile} -Bx${axistick[0]}+l\"${labels[0]}\" \
-By${axistick[1]}+l\"${labels[1]}\" -JX${layouts[0]}/${pic_height}i \
-X${image_origin[0]}i -Y${image_origin[1]}i -K -P ${from_open} $psfile"
else
gmt grdimage ${data} -R${range} -C${cptfile} -Bx${axistick[0]}+l"${labels[0]}" \
-By${axistick[1]}+l"${labels[1]}" -JX${layouts[0]}/${pic_height}i \
-X${image_origin[0]}i -Y${image_origin[1]}i -K -P ${from_open} $psfile
fi
fi
# 画多边形
if [[ $polyfile != 'null' ]]; then
RunOrEcho ${run_command} "gmt psxy ${polyfile} -JX${layouts[0]}/${pic_height}i -W0.25p,black,- -R${range} -L -K -O >> $psfile"
fi
#添加号码
if [[ ${annotation} != 'null' ]]; then
if [[ ${run_command} == 0 ]]; then
echo "gmt pstext -R${range} -JX${layouts[0]}/${pic_height}i -K -O <<- EOF >> $psfile"
echo ${annotation}
echo "EOF"
else
gmt pstext -R${range} -JX${layouts[0]}/${pic_height}i -K -O <<- EOF >> $psfile
${annotation}
EOF
fi
fi
#-C${cptfile}+Uk 使用km色标单位除1000
#如果unit等于km则在cptfile后面添加+Uk
if [[ ${plot_colorbar} == 1 ]]; then
if [[ ${unit} == 'km+Uk' ]]; then
RunOrEcho ${run_command} "gmt psscale -Dx${layouts[2]}+w${layouts[3]} -C${cptfile}+Uk -Bx${bartick} -By+l${unit} \
-X${bar_origin[0]}i -Y${bar_origin[1]}i -O -K >> $psfile"
else
RunOrEcho ${run_command} "gmt psscale -Dx${layouts[2]}+w${layouts[3]} -C${cptfile} -Bx${bartick} -By+l${unit} \
-X${bar_origin[0]}i -Y${bar_origin[1]}i -O -K >> $psfile"
fi
fi
if [[ ${keep_open} == 0 ]]; then
# 这一句什么都不做 唯一的作用就是关闭ps文件
if [[ ${run_command} == 0 ]]; then
echo "gmt pstext -R${range} -JX -O <<- EOF >> $psfile"
echo "EOF"
else
gmt pstext -R${range} -JX -O <<- EOF >> $psfile
EOF
fi
# 输出 eps 和 png 文件
RunOrEcho ${run_command} "gmt psconvert $psfile -A -TEG -E300"
# 删除临时文件 使用linux终端rm命令
RunOrEcho ${run_command} "rm $cptfile $psfile gmt.history gmt.conf"
if [[ $plotgrad == 1 ]]; then
RunOrEcho ${run_command} "rm $gradfile"
fi
# 在终端显示图像 此命令需要imgcat.sh脚本和iTerm终端
#imgcat $jpgfile
# 打开图片文件 此命令使用MacOS终端open命令
if [[ ${run_command} == 1 ]]; then
open $jpgfile
fi
fi
fi