------------------------------
List of the following files:
------------------------------
jcope2.run.sh
jcope2_sst.grad.f90
jcope2_ysec.f90
makefile
plot.jcope2_sst.grad.sh
gmtpar.sh
note.sh
------------------------------
Machine info
------------------------------
aofd165.bio.mie-u.ac.jp
/work2/am/Kuroshio.Kii.Hantou/JCOPE2.SST
Sat Sep 23 16:12:21 JST 2017
======================
jcope2.run.sh
======================
#!/bin/sh
exe1=jcope2_sst.grad.run.sh
exe2=jcope2_ysec.run.sh
plt=plot.jcope2_sst.grad.sh
exefiles="$exe1 $exe2 $plt"
for ef in $exefiles; do
if [ ! -f $ef ]; then
echo
echo ERROR in $0: No such file, $ef
echo
exit 1
fi
done
# time info
syyyy=2000
eyyyy=2009
smm=06
sdd=01
days=30 #30
lonout=136.1
yyyy=$syyyy
while [ $yyyy -le $eyyyy ]; do #yyyy
d=0
while [ $d -lt $days ]; do
yyyymmdd=$(date -d"${yyyy}${smm}${sdd} $d day" '+%Y%m%d')
$exe1 $yyyymmdd
$exe2 $yyyymmdd $lonout
$plt $yyyymmdd $lonout
d=$(expr $d + 1)
done #d
yyyy=$(expr $yyyy + 1)
done #yyyy
exit 0
----------------------
End of jcope2.run.sh
----------------------
======================
jcope2_sst.grad.f90
======================
program jcope2_sst
real,allocatable ::tj(:,:,:),uj(:,:,:),vj(:,:,:),el(:,:) !JCOPE2 data
real, allocatable :: lon(:),lat(:)
real, allocatable :: isst1(:,:),u(:,:),v(:,:)
real, allocatable :: sstgrad(:,:) !gradient
real,parameter::a=6378.137e3, pi=3.14159265359, dtr=pi/180.0, &
r2d=1.0/dtr
real,parameter::m_to_100km=1.e5, m2cm=100.0
! other variables
character(256) :: ifile1, ifile2, ifile3, ifile4,ofile
character(80) :: cxfcst ! String of xfcst
integer, parameter :: iunit = 11, ounit = 21
integer yyyy,mm,dd
real ,parameter::fillvalue=-999.0
character(4) :: fildsc
integer iyf,imf,idf,mall
namelist /para/ifile1, ifile2, ifile3, ifile4, ofile, &
yyyy, mm, dd, im, jm, km, tbias, &
xlons, xlons_offset_denom, ylats, ylats_offset_denom, deltalat, &
deltalon, iouts,ioute,jouts,joute
read(*,para)
allocate(tj(im,jm,km),uj(im,jm,km),vj(im,jm,km),el(im,jm))
allocate(lon(im),lat(jm),isst1(im,jm),u(im,jm),v(im,jm))
allocate(sstgrad(im,jm))
startlon =xlons-1.0/xlons_offset_denom
startlat =ylats-1.0/ylats_offset_denom
! read JCOPE2 Temp data
tj(:,:,:)=0.0
uj(:,:,:)=0.0
vj(:,:,:)=0.0
print '(A,A)','Reading ',trim(ifile1)
open(iunit, file=ifile1,form='unformatted',action="read")
read(iunit)&
& fildsc,iyf,imf,idf &
& ,(((tj(i,j,k),i=1,im),j=1,jm),k=1,km-1),mall
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,iyf,imf,idf,mall = ',fildsc,iyf,imf,idf,mall
close(iunit)
print '(A)','Done.'
print '(A,A)','Reading ',trim(ifile2)
open(iunit, file=ifile2,form='unformatted',action="read")
read(iunit)&
& fildsc,iyf,imf,idf &
& ,(((uj(i,j,k),i=1,im),j=1,jm),k=1,km-1),mall
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,iyf,imf,idf,mall = ',fildsc,iyf,imf,idf,mall
close(iunit)
print '(A)','Done.'
print '(A,A)','Reading ',trim(ifile3)
open(iunit, file=ifile3,form='unformatted',action="read")
read(iunit)&
& fildsc,iyf,imf,idf &
& ,(((vj(i,j,k),i=1,im),j=1,jm),k=1,km-1),mall
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,iyf,imf,idf,mall = ',fildsc,iyf,imf,idf,mall
close(iunit)
print '(A)','Done.'
print '(A,A)','Reading ',trim(ifile4)
open(iunit, file=ifile4,form='unformatted',action="read")
read(iunit)&
& fildsc,iyf,imf,idf &
& ,((el(i,j),i=1,im),j=1,jm),mall
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,iyf,imf,idf,mall = ',fildsc,iyf,imf,idf,mall
close(iunit)
isst1(:,:)=tj(:,:,3)+tbias
where(isst1 == tbias) isst1 = fillvalue
u(:,:)=uj(:,:,3) !*m2cm
v(:,:)=vj(:,:,3) !*m2cm
do i=1,im
lon(i)=startlon+float(i-1)*deltalon
enddo !j
do j=1,jm
lat(j)=startlat+float(j-1)*deltalat
enddo !j
sstgrad=fillvalue
dy=a*dtr*deltalat*2.0
do j=jouts,joute !2,jm-1
dx=a*cos(lat(j)*dtr)*deltalon*2.0
do i=iouts,ioute !2,im-1
sstgrad_x=(isst1(i+1,j) -isst1(i-1,j ))/dx
sstgrad_y=(isst1(i ,j+1)-isst1(i ,j-1))/dy
if(isst1(i,j)/=fillvalue .and. isst1(i-1,j)/=fillvalue .and. &
isst1(i+1,j)/=fillvalue .and. isst1(i,j-1)/=fillvalue .and. &
isst1(i,j+1)/=fillvalue )then
sstgrad(i,j)=sqrt(sstgrad_x**2 + sstgrad_y**2) * m_to_100km
else
sstgrad(i,j)=fillvalue
endif
end do !i
end do !j
! write data
open(ounit, file=ofile) !, form='unformatted', access='sequential')
do j=jouts,joute
do i=iouts,ioute
ang=atan2(v(i,j),u(i,j))*r2d
vmg=sqrt(u(i,j)**2+v(i,j)**2)
write(ounit,'(f9.4,f8.4,1x,f7.2,e12.4,2x,2f7.2,1x,f7.3)')&
lon(i),lat(j),isst1(i,j),sstgrad(i,j),ang,vmg,el(i,j)
enddo !i
enddo !j
print *
print '(A,A)','Output: ',trim(ofile)
print *
end program jcope2_sst
----------------------
End of jcope2_sst.grad.f90
----------------------
======================
jcope2_ysec.f90
======================
program jcope2_ysec
real,allocatable ::tj(:,:,:),uj(:,:,:),vj(:,:,:),sj(:,:,:) !JCOPE2 data
real, allocatable :: lon(:),lat(:)
real, allocatable :: t(:,:),u(:,:),v(:,:),s(:,:)
real, allocatable :: z(:,:,:),zz(:,:,:),dz(:,:,:)
real, allocatable :: sstgrad(:,:) !gradient
! VERTICAL INTERPOLATION
real,allocatable :: h0(:)
real,allocatable ::tlev(:,:,:),ulev(:,:,:),vlev(:,:,:),slev(:,:,:)
real,parameter::a=6378.137e3, pi=3.14159265359, dtr=pi/180.0, &
r2d=1.0/dtr
real,parameter::m_to_100km=1.e5, m2cm=100.0
integer,parameter::klev=41
real lonout,latouts, latoutn
! other variables
character(256) :: ifile1, ifile2, ifile3, ifile4, ifile5, ofile
character(80) :: cxfcst ! String of xfcst
integer, parameter :: iunit = 11, ounit = 21
integer yyyy,mm,dd
real,parameter::fillvalue=-999.0, valmis=fillvalue
real,parameter :: tbias=10.0, sbias=35.0
character(4) :: fildsc
integer iyf,imf,idf,mall
namelist /para/ifile1, ifile2, ifile3, ifile4, ifile5, ofile, &
yyyy, mm, dd, im, jm, km, lonout, latouts, latoutn, &
xlons, xlons_offset_denom, ylats, ylats_offset_denom, deltalat, &
deltalon
read(*,para)
allocate(tj(im,jm,km),uj(im,jm,km),vj(im,jm,km),sj(im,jm,km))
allocate(z(im,jm,km),zz(im,jm,km),dz(im,jm,km))
allocate(lon(im),lat(jm),t(jm,km),u(jm,km),v(jm,km),s(jm,km))
allocate(h0(klev))
zinc=5.0
do k=1,klev
h0(k)=0.0-zinc*float(k-1)
end do !k
allocate(ulev(im,jm,klev),vlev(im,jm,klev) &
&,tlev(im,jm,klev),slev(im,jm,klev))
startlon =xlons-1.0/xlons_offset_denom
startlat =ylats-1.0/ylats_offset_denom
! read JCOPE2 Temp data
tj(:,:,:)=0.0
uj(:,:,:)=0.0
vj(:,:,:)=0.0
sj(:,:,:)=0.0
print '(A,A)','Reading ',trim(ifile1)
open(iunit, file=ifile1,form='unformatted',action="read")
read(iunit)&
& fildsc,iyf,imf,idf &
& ,(((tj(i,j,k),i=1,im),j=1,jm),k=1,km-1),mall
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,iyf,imf,idf,mall = ',fildsc,iyf,imf,idf,mall
close(iunit)
print '(A)','Done.'
print '(A,A)','Reading ',trim(ifile2)
open(iunit, file=ifile2,form='unformatted',action="read")
read(iunit)&
& fildsc,iyf,imf,idf &
& ,(((uj(i,j,k),i=1,im),j=1,jm),k=1,km-1),mall
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,iyf,imf,idf,mall = ',fildsc,iyf,imf,idf,mall
close(iunit)
print '(A)','Done.'
print '(A,A)','Reading ',trim(ifile3)
open(iunit, file=ifile3,form='unformatted',action="read")
read(iunit)&
& fildsc,iyf,imf,idf &
& ,(((vj(i,j,k),i=1,im),j=1,jm),k=1,km-1),mall
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,iyf,imf,idf,mall = ',fildsc,iyf,imf,idf,mall
close(iunit)
print '(A)','Done.'
print '(A,A)','Reading ',trim(ifile4)
open(iunit, file=ifile4,form='unformatted',action="read")
read(iunit)&
& fildsc,iyf,imf,idf &
& ,(((sj(i,j,k),i=1,im),j=1,jm),k=1,km-1),mall
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,iyf,imf,idf,mall = ',fildsc,iyf,imf,idf,mall
close(iunit)
open(iunit,file=ifile5,form='unformatted',action="read")
read(iunit) fildsc,Z,ZZ,DZ,ichflg
print '(A,A4,1x,i4.4,1x,i2.2,1x,i2.2,1x,i5)', &
'fildsc,ichflg = ',fildsc,ichflg
close(iunit)
do i=1,im
lon(i)=startlon+float(i-1)*deltalon
enddo !j
do j=1,jm
lat(j)=startlat+float(j-1)*deltalat
enddo !j
idx=int((lonout-startlon)/deltalon)+2
print '(A,f9.3,A,i4,A,f9.3)','lonout= ',lonout, ' idx= ',idx,' lon(idx)=',lon(idx)
jdxs=int((latouts-startlat)/deltalat)+1
jdxn=int((latoutn-startlat)/deltalat)+1
print '(A,f9.3,A,i4,A,f9.3)','latouts=',latouts, ' jdxs=',jdxs,' lat(jdxs)=',lat(jdxs)
print '(A,f9.3,A,i4,A,f9.3)','latoutn=',latoutn, ' jdxn=',jdxn,' lat(jdxn)=',lat(jdxn)
sstgrad=fillvalue
dy=a*dtr*deltalat*2.0
do j=jouts,joute !2,jm-1
dx=a*cos(lat(j)*dtr)*deltalon*2.0
end do !j
do k=1,klev
call gsigtoz3 &
& ( zz,uj,ulev,h0(k),im,jm,km,klev,k,valmis )
call gsigtoz3 &
& ( zz,vj,vlev,h0(k),im,jm,km,klev,k,valmis )
call gsigtoz3 &
& ( zz,tj,tlev,h0(k),im,jm,km,klev,k,valmis )
call gsigtoz3 &
& ( zz,sj,slev,h0(k),im,jm,km,klev,k,valmis )
enddo
do i=1,im
do j=1,jm
do k=1,klev
if( abs(z(i,j,km)).le.1. ) then
ulev(i,j,k) = valmis
vlev(i,j,k) = valmis
tlev(i,j,k) = valmis
slev(i,j,k) = valmis
else
tlev(i,j,k) = tlev(i,j,k) +tbias
if (tlev(i,j,k).lt.-20.0) tlev(i,j,k)=valmis
slev(i,j,k) = slev(i,j,k) +sbias
if (slev(i,j,k).lt.-20.0) slev(i,j,k)=valmis
end if
enddo
enddo
enddo
! write data
open(ounit, file=ofile) !, form='unformatted', access='sequential')
write(ounit,'(A,f10.5)')'# lonout= ', lon(idx)
write(ounit,'(A,f10.5)')'# latouts= ',lat(jdxs)
write(ounit,'(A,f10.5)')'# latoutn= ',lat(jdxn)
write(ounit,'(A)')'# lat(j),y,h0(k),tlev(idx,j,k),slev(idx,j,k),ulev(idx,j,k),vlev(idx,j,k)'
y0=dy*float(jdxs)/1000.0
do j=jdxs,jdxn
y=dy*float(j)/1000.0-y0
do k=1,klev
write(ounit,'(f9.3,f7.1,f10.2,f9.2,f10.3,1x,2f8.2)')&
lat(j),y,h0(k),tlev(idx,j,k),slev(idx,j,k),ulev(idx,j,k),vlev(idx,j,k)
end do !k
end do !j
print *
print '(A,A)','Output: ',trim(ofile)
print *
end program jcope2_ysec
SUBROUTINE GSIGTOZ3(ZZ,T,TLEV,ZLEV,IM,JM,KM,klev,kk,valmis)
!C-------------------------------------------------------------------
!C THIS ROUTINE LINERLY INTERPOLATES TLEV AT THE LEVEL,
!C ZLEV, FROM T LOCATED ON SIGMA LEVELS, ZZ.
!C ZLEV AND ZZ ARE NEGATIVE QUANTITIES.
!C A SEARCH IS MADE TO FIND ZZ(I,J,K) AND ZZ(I,J,K+1) WHICH
!C BRACKETS ZLEV; THEN THE INTERPOLATION IS MADE.
!C IN THE REGION < 0 and > ZZ(I,J,1) AND, IN THE REGION,
!C < ZZ(I,J,KM-1) and > -1, DATA IS EXTRAPOLATED.
!C NOTE THAT A NEW MASK ,FSM, IS CREATED.
!C
!C THE VALUES OF H SUPPLIED TO THIS SUBROUTINE SHOULD BE
!C APPROPRIATE TO THE VARIABLE T. FOR EXAMPLE, IF T IS THE
!C X-COMPONENT OF VELOCITY, H SHOULD THE AVERAGE OF DEPTH
!C AT I AND I-1.
!C-------------------------------------------------------------------
DIMENSION ZZ(IM,JM,KM),T(IM,JM,KM),TLEV(IM,JM,KLEV)
do j=1,jm
do i=1,im
tlev(i,j,kk) = valmis
enddo
enddo
DO 200 J=1,JM
DO 200 I=1,IM
IF(ZLEV.GT.ZZ(I,J,1) .and. t(i,j,1).gt.valmis+1. .and. &
t(i,j,2).gt.valmis+1. ) THEN
K=1
TLEV(I,J,kk)=T(I,J,K)+(ZLEV-ZZ(I,J,K)) &
*(T(I,J,K+1)-T(I,J,K))/(ZZ(I,J,K+1)-ZZ(I,J,K))
GO TO 200
ENDIF
K=1
100 CONTINUE
IF(ZLEV.LE.ZZ(I,J,K).AND.ZLEV.GE.ZZ(I,J,K+1) .and. &
t(i,j,k).gt.valmis+1. .and. t(i,j,k+1).gt.valmis+1. ) THEN
TLEV(I,J,kk)=T(I,J,K)+(ZLEV-ZZ(I,J,K)) &
*(T(I,J,K+1)-T(I,J,K))/(ZZ(I,J,K+1)-ZZ(I,J,K))
GO TO 200
ELSE
K=K+1
IF(K.LT.(KM-1)) GO TO 100
ENDIF
IF(ZLEV.LE.ZZ(I,J,K).AND.ZLEV.GE.ZZ(I,J,KM) .and. &
t(i,j,k).gt.valmis+1. .and. t(i,j,k-1).gt.valmis+1. ) THEN
TLEV(I,J,kk)=T(I,J,K-1)+(ZLEV-ZZ(I,J,K-1)) &
*(T(I,J,K)-T(I,J,K-1))/(ZZ(I,J,K)-ZZ(I,J,K-1))
ENDIF
200 CONTINUE
RETURN
END
----------------------
End of jcope2_ysec.f90
----------------------
======================
makefile
======================
#
#
#j Fortran77とFortran90/95の混成プログラムのコンパイル用makefile
#
#j 作成 Sun Feb 15 15:16:02 2009
#
#j 注:ごく簡単な動作チェックしか行っていない
#j サフィックスルールがこれで正しいかどうかチェックする必要あり
#
#j 参考文献
#j makefile の書き方
#j http://www.eis.osakafu-u.ac.jp/~yabu/soft/makefile.html
#j Makefile小技と新バージョンのDelFEM(FEMとUIの日記@New York)
#j http://d.hatena.ne.jp/etopirika5/20091207/1260207955
#
#-----------------------------------------
#j マクロの定義
#-----------------------------------------
#j コンパイラの指定
FC=ifort
#FC=gfortran
#FC=g77
#FC=f77
#j OBJDIR (オブジェクトファイルをおくディレクトリ)
OBJDIR=./obj
#-----------------------------------------
#j コンパイルオプション(ifort)
#-----------------------------------------
#j デバッグ用オプション(はじめて実行するときには必ずこのオプションをつけてコンパイルする)
#FDFLAGS= -CB -traceback -fpe0
#j コンパイルオプション
FFLAGS=-module ${OBJDIR}
#j -module : モジュールファイル(.mod)を置く場所を指定する
#j 最適化用(他にも何種類かある. ifort -helpかマニュアルを見て調べる)
#FFLAGS = -O3 -module ${OBJDIR}
#j 倍精度
#FFLAGS = -r8 -module ${OBJDIR}
#j 入力データのバイト・スワップ(大型計算機のバイナリデータを読む時などにつかう)
FFLAGS = -O3 -convert big_endian -module ${OBJDIR}
#-----------------------------------------
#j リンク用のオプション
#-----------------------------------------
LDFLAGS=-module ${OBJDIR}
#j ターゲット名(最終的に作りたい実行ファイル名)
TARGET1=jcope2_ysec.exe
TARGET2=jcope2_sst.grad.exe
TARGETS=$(TARGET1) $(TARGET2)
#j MOD (Fortran90以降のモジュールファイル)
MOD=
#j SRC7 (Fortran77のソースファイル)
SRC7=
#j SRC9 (Fortran90のソースファイル)
SRC9=jcope2_ysec.f90
#j OBJ1 (ターゲットを作るのに必要なオブジェクトファイル名)
OBJM=$(MOD:.f90=.o)
OBJ7=$(SRC7:.f=.o)
OBJ9=$(SRC9:.f90=.o)
OBJTMP=$(OBJM) $(OBJ7) $(OBJ9)
OBJ1=$(OBJTMP:%=${OBJDIR}/%)
MOD2=
#j SRC7 (Fortran77のソースファイル)
SRC72=
#j SRC9 (Fortran90のソースファイル)
SRC92=jcope2_sst.grad.f90
#j OBJ (ターゲットを作るのに必要なオブジェクトファイル名)
OBJM2=$(MOD2:.f90=.o)
OBJ72=$(SRC72:.f=.o)
OBJ92=$(SRC92:.f90=.o)
OBJTMP2=$(OBJM2) $(OBJ72) $(OBJ92)
OBJ2=$(OBJTMP2:%=${OBJDIR}/%)
OBJ=$(OBJ1) $(OBJ2)
#-----------------------------------------
#j ここからコンパイルのルールの記述
#-----------------------------------------
all: mkobjd $(TARGETS)
$(TARGET1): $(OBJ1)
$(FC) -o $@ $(OBJ1) ${LDFLAGS}
$(TARGET2): $(OBJ2)
$(FC) -o $@ $(OBJ2) ${LDFLAGS}
#j 暗黙のサフィックスルールを無効にする
.SUFFIXES:
#j サフィックスの追加
.SUFFIXES: .o .f .f90
#j f77のソースファイルのコンパイル(.f用のサフィックスルール)
${OBJDIR}/%.o : %.f
$(FC) -c ${FDFLAGS} ${FFLAGS} -c -o $@ $<
#j f90のソースファイルのコンパイル(.f90用のサフィックスルール)
${OBJDIR}/%.o : %.f90
$(FC) -c ${FDFLAGS} ${FFLAGS} -c -o $@ $<
#j オブジェクトファイルを削除
clean:
rm -rf $(OBJ) ${OBJDIR}/*.mod
#j オブジェクトファイルと実行ファイルを削除
distclean:
rm -rf $(OBJ) ${OBJDIR}/*.mod $(TARGETS) \
rm -rf ${OBJDIR}
#j オブジェクトファイルを置くディレクトリを作成
mkobjd:
if [ ! -d ${OBJDIR} ]; then \
mkdir -p ${OBJDIR} ; \
fi
run:
cd .. ; \
./run.sh 2>&1 > log.txt ; \
cd src
#j tarファイルを作る
tar:
tar cvf $(TARGET1).tar ./*.f90 ./*.f ./*.c ./*.h ./*.txt makefile
----------------------
End of makefile
----------------------
======================
plot.jcope2_sst.grad.sh
======================
#!/bin/bash
# Description:
#
# Author: am
#
# Host: aofd165.bio.mie-u.ac.jp
# Directory: /work2/am/Kuroshio.Kii.Hantou/JCOPE2.SST
#
# Revision history:
# This file is created by /usr/local/mybin/ngmt.sh at 10:32 on 09-23-2017.
. ./gmtpar.sh
gmtset HEADER_FONT_SIZE 16
gmtset LABEL_FONT_SIZE 14
gmtset ANOT_FONT_SIZE 14
echo "Bash script $0 starts."
range=134/140/31/35
size=Q137/4
xanot=a2f1
yanot=a2f1
anot=${xanot}/${yanot}WsNe
reso=0.05/0.05
tension=0.8
yyyymmdd=$1
lonout=$2
indir=output
in1=${indir}/sst.grad_${yyyymmdd}.txt
in2=${indir}/ysec.lon${lonout}_${yyyymmdd}.txt
infiles="$in1 $in2"
for input in $infiles; do
if [ ! -f $input ]; then
echo Error in $0 : No such file, $input
exit 1
fi
done
figdir="fig"
if [ ! -d ${figdir} ];then
mkdir -p $figdir
fi
out=${figdir}/$(basename $in1 .txt).ps
grd=$(basename $in1 .txt).grd
cpt=$(basename $in1 .txt).cpt
awk '{if ($3>-2.0) print $1,$2,$3}' $in1|minmax
makecpt -Cgray -I -T2/10/2 -Z > $cpt
awk '{if ($4>0.0) print $1,$2,$4}' $in1|\
surface -R$range -I$reso -T$tension -G${grd}
landmask=landmask.grd
maskedgrd=masked.grd
grdlandmask -R$range -Dl -I$reso -N1/NaN -G$landmask #-V
grdmath $grd $landmask OR = $maskedgrd
grdimage $maskedgrd -R$range -J$size -C$cpt -K -Y7 -P > $out
psscale -D1.5/-0.2/2.5/0.075h -C$cpt -B2:"":/:K/100km: -O -K >> $out
# EL
#awk '{if ($7>-100.0) print $1,$2,$7}' $in1|\
#surface -R$range -I$reso -T$tension -G${grd}
#grdcontour $grd -R$range -J$size -W2${dash} -C0.2 -A0.2f10 -O -K >> $out
# ISOTACH
awk '{if ($6>-100.0) print $1,$2,$6}' $in1|\
surface -R$range -I$reso -T$tension -G${grd}
grdcontour $grd -R$range -J$size -W4${dash} -L1/100 -C0.5 -A0.5f10 -O -K >> $out
# Vector Size
#VSIZE=1.5 #0.1 # ms-1
#VSCALE=0.2 # inch
#VFACTOR=`echo "scale=6; ${VSCALE}/${VSIZE}" | bc`
#
#awk -v SC="$VFACTOR" '{if($6>=1.0 && NR%2 ) \
#{printf "%12.6f %12.6f %10.3f %10.4f\n", $1, $2, \
#$5, $6*SC} }' $in | psxy -R -J${size} \
# -Sv0.01/0.1/0.03n0.2 -G0 \
# -O -K >> $out
awk '{if ($3>-2.0) print $1,$2,$3}' $in1|\
surface -R$range -I$reso -T$tension -G${grd}
grdcontour $grd -R$range -J$size -W2 -C.5 -A1f10 -O -K >> $out
pscoast -R$range -J$size -B$anot -W3 -G255 -Di -O -K >>$out
lonr=$(awk '{if($2=="lonout=")print $3}' $in2)
lats=$(awk '{if($2=="latouts=")print $3}' $in2)
latn=$(awk '{if($2=="latoutn=")print $3}' $in2)
echo "lonr lats latn = $lonr $lats $latn"
psxy <<EOF -R$range -J$size -W8 -O -K >>$out
$lonr $lats
$lonr $latn
EOF
pstext <<EOF -R$range -J$size -D-0.05/0 -Wwhite -O -K >>$out
$lonr $lats 14 0 0 MR A
$lonr $latn 14 0 0 MR B
EOF
#
# Y SECTION
#
reso2=0.05/5
grd2=$(basename $in2 .txt).grd
range2=31.95/34.05/0/200
range2out=32/34/0/200
size2=2/-2.66
cpt2=$(basename $in2 .txt).cpt
makecpt -Cgray -I -T5/30/1 -Z > $cpt2
awk '{if ($4>-2.0) print $1,-$3,$4}' $in2|minmax
awk '{if ($4>-2.0) print $1,-$3,$4}' $in2|\
surface -R$range2 -I$reso2 -T$tension -G${grd2}
grdimage $grd2 -R$range2out -JX$size2 -C$cpt2 -O -K -X4.4 >> $out
grdcontour $grd2 -R$range2out -JX$size2 -C1 -W2/${white} -O -K >>$out
psscale -D1/-0.2/1.5/0.075h -C$cpt2 -B5f5:"":/:${deg}C: -O -K >> $out
awk '{if ($6>-2.0) print $1,-$3,$6}' $in2|minmax
awk '{if ($6>-2.0) print $1,-$3,$6}' $in2|\
surface -R$range2 -I$reso2 -T$tension -G${grd2}
grdcontour $grd2 -R$range2out -JX$size2 -W10/${white} -L0/100 -C0.5 -O -K >> $out
grdcontour $grd2 -R$range2out -JX$size2 -W3 -L0/100 -C0.5 -A0.5f10 -G1/2 -O -K >> $out
grdcontour $grd2 -R$range2out -JX$size2 -W10/${white} -L-100/-0.5 -C0.5 -O -K >> $out
grdcontour $grd2 -R$range2out -JX$size2 -W3${dashed} -L-100/-0.5 -C0.5 -G1/2 -A0.5f10 -O -K >> $out
psbasemap -R$range2out -JX$size2 -Ba1f0.5:"Latitude":/a50:Depth${sp}[m]:wsNE -O -K >>$out
#ulegmag=$VSIZE
#uscl=${VSCALE}
#ufac=$(echo "scale=5; $uscl/$ulegmag" | bc)
#
range2=0/1/0/1
size2=8/1
#legend=legend.txt
#cat <<EOF >$legend
#0.6 0.35 0.0 ${ulegmag}
#EOF
#awk -v f=$ufac '{if($1!="#") \
#printf "%12.6f %12.6f %12.6f %12.6f\n", \
#$1,$2,$3,$4*f}' $legend| \
#psxy -R$range2 -JX$size2 -Sv0.01/0.1/0.05n0.2 -G0 \
# -K -O -Y3 >> $out
#echo "0.65 0.35 16 0 0 ML ${ulegmag} m/s" |pstext -R -JX -O -K >> $out
pstext <<EOF -R$range2 -JX$size2 -O -K -X-4.2 -Y2.8 >>$out
0.52 0.3 14 0 0 MC A
0.77 0.3 14 0 0 MC B
EOF
pstext <<EOF -R$range2 -JX$size2 -O -K -X-1 >>$out
0.1 0.53 18 0 0 MC (a)
0.64 0.53 18 0 0 MC (b)
EOF
yyyy=${yyyymmdd:0:4}
m=${yyyymmdd:4:2}
dd=${yyyymmdd:6:7}
if [ $m -eq 1 ]; then mmm="Jan"; fi
if [ $m -eq 2 ]; then mmm="Feb"; fi
if [ $m -eq 3 ]; then mmm="Mar"; fi
if [ $m -eq 4 ]; then mmm="Apr"; fi
if [ $m -eq 5 ]; then mmm="May"; fi
if [ $m -eq 6 ]; then mmm="Jun"; fi
if [ $m -eq 7 ]; then mmm="Jul"; fi
if [ $m -eq 8 ]; then mmm="Aug"; fi
if [ $m -eq 9 ]; then mmm="Sep"; fi
if [ $m -eq 10 ]; then mmm="Oct"; fi
if [ $m -eq 11 ]; then mmm="Nov"; fi
if [ $m -eq 12 ]; then mmm="Dec"; fi
echo ${dd}${mmm}${yyyy}
pstext <<EOF -R$range2 -JX$size2 -O -K >>$out
0.5 0.8 18 0 0 MC ${dd}${mmm}${yyyy}
EOF
jpg=$(basename $out .ps).jpg
xoffset=
yoffset=0.7
in="$infiles"
comment="Output : ${out}" #${jpg}"
export LANG=C
. ./note.sh
echo Convert $out to ${figdir}/$jpg
convert $out ${figdir}/$jpg
echo
rm -rf $grd $maskedgrd $cpt
rm -f $grd2 $cpt2
echo
echo Input : $in1
echo Input : $in2
echo Output :$out # $jpg #$out
echo
echo "Done $0"
----------------------
End of plot.jcope2_sst.grad.sh
----------------------
======================
gmtpar.sh
======================
#
# GMTのパラメータの設定 (version 4以上対象)
#
# 長さの単位はインチとする
gmtset MEASURE_UNIT INCH
gmtset LABEL_FONT_SIZE 18
gmtset HEADER_FONT_SIZE 20
gmtset ANOT_FONT_SIZE 18
gmtset TICK_PEN 4
# 地図の縦横軸に縞々を使わない
gmtset BASEMAP_TYPE PLAIN
# 紙のサイズはA4
gmtset PAPER_MEDIA A4
# 地図の°の記号
gmtset DEGREE_SYMBOL = degree # Available for ver. 4 or later
# 空白文字の設定
sp="\040" # White space
# =の記号
eq="\075" # equal
# 温度の°の記号
deg="\260" #deg="\312" # degree symbol
# 色のRGB値を設定 (線や点に色をつけるときRGB値を直接指定するより分かりやすい)
red="255/0/0"
orange="255/165/0"
pink="255/192/203"
green="0/255/0"
blue="0/0/255"
midnightblue="25/25/112"
yellow="255/255/0"
gold="255/215/0"
purple="160/32/240"
magenta="255/0/255"
white="255/255/255"
# 線種を指定
dash="t15_5:15"
dot="t5_5:5"
dotdash="t12_5_5_5:5"
#その他の線種の例
# -W5t20_5:5
# -W5t15_5:5
# -W5t20_5_5_5_5_5:5
# -W5t30_5_5_5:5
# -W5t20_5_10_5:5
# 数字の出力フォーマット(project, grd2xzy等で使う)
gmtset D_FORMAT %lg #デフォルト(規定値)
# 桁数を指定(例:123.45678)
#gmtset D_FORMAT %12.6f
# 状況に応じて
#gmtset D_FORMAT %12.5g
----------------------
End of gmtpar.sh
----------------------
======================
note.sh
======================
# Print time, current working directory and output filename
export LANG=C
currentdir=`pwd`
if [ ${#currentdir} -gt 90 ]; then
curdir1=${currentdir:1:90}
curdir2=${currentdir:91}
else
curdir1=$currentdir
curdir2="\ "
fi
now=`date`
host=`hostname`
#comment=" "
time1=$(ls -l $in | awk '{print $6, $7, $8}')
pstext -JX6/1.2 -R0/1/0/1.2 -N -X${xoffset:-0} -Y${yoffset:-0} << EOF -O >> $out
0 1.1 9 0 1 LM $0 $@
0 0.95 9 0 1 LM ${now}
0 0.80 9 0 1 LM ${host}
0 0.65 9 0 1 LM ${curdir1}
0 0.50 9 0 1 LM ${curdir2}
0 0.35 9 0 1 LM Input: ${in} (${time1})
0 0.1 9 0 1 LM ${comment:-""}
EOF
# Output: ${out}
----------------------
End of note.sh
----------------------