各海域における海面高度の上昇傾向
今回はCMEMSが配布する海面高度データを使用する。
海面温度データと同じく、３６０ｘ１８０の1度グリッドデータとしている。
以下の図は、季節変動を除去したのちに、線形トレンドを計算してプロットしている。
なお、1993年～2000年までの上昇トレンドは1度グリッドでは2.66ｍｍ/年となった。
2013年以降の上昇トレンドを計算してみよう。
似たような計算結果はここでもみられる。ただしグリッドサイズが異なるので結果は若干違う。
計算の手順
データはここ！
１）1993年1月～2023年5月までの各月における平均海面水位（実際には、1993年～2012年の平均場からのズレとなる海面偏差の平均）を計算する。→平均海面水位の時系列の取得
２）平均海面水位を従属変数、時間を説明変数として線形一次の回帰式を求める。
Ｙ（海面水位）＝a x X(時間, xvar)+b
３）時間xvarと、求めた傾きaと切片b、で予測される1993年～2023年までの海面水位を求める。
４）元の時系列に季節変動成分が含まれている。各月毎に、１）の時系列と３）の時系列の差の平均を求める。
サンプルコードでは変数monmeanとしている。
５）１）の時系列から４）を引くことで、季節変動シグナルを除去した海面水位時系列を得る。
６）５）の時系列に対して、任意の期間における線形一次の回帰式を求める。
得られた傾きはひと月での上昇率となるので、１２を掛けて1年での上昇率に変換する。
サンプルコード
c234567--------------------------------
character*23 ifl
integer inum,inum2,smo,emo
real dummy(360,180),msla
real sla3d(360,180,365),xvar(365)
real sla1d(360),sla1d2(365)
real sla1df(365),sla1df2(365)
real sla1d9303(120)
real sla1d1323(125)
real slope1,offset1
real slope2,offset2
real slope3,offset3
integer month(365),year(365)
real monmean(12)
!--------------------------------------
! read SLA dataset; build SLA 3d(lon,lat,time)
!--------------------------------------
inum=1
smo=1;emo=12;
do iyr=1993,2023
if(iyr.eq.2023)then
emo=5
endif
do imo=smo,emo
if(imo.ge.10)then
write(ifl,'("./CMEMS/",i4,"_m",i2,"sla.bin")') iyr,imo
else
write(ifl,'("./CMEMS/",i4,"_m0",i1,"sla.bin")') iyr,imo
endif
write(*,200) ifl
200 format(a23) !--make file's name
open(10,file=ifl,form='unformatted',
& access="direct",status="old",recl=360*180*4)
read(10,rec=1) dummy
do j=1,180
do i=1,360
sla3d(i,j,inum)=dummy(i,j)
enddo
enddo
year(inum)=iyr ! we will use it later!
month(inum)=imo ! we will use it later!
xvar(inum)=inum
! write(*,*) inum,year(inum),month(inum),xvar(inum)
inum=inum+1
enddo
enddo
inum=inum-1
write(*,*) "inum=",inum
!--------------------------------------
! cal mean value of timeseries
!--------------------------------------
do k=1,inum
sla1d(k)=0;
hoge=0
inum2=0;
do j=1,180 !summation!
do i=1,360
if(sla3d(i,j,k).gt.-900)then
!sla1d(k)=sla1d(k)+sla3d(i,j,k)
hoge=hoge+sla3d(i,j,k)
inum2=inum2+1
endif
enddo
enddo
sla1d(k)=hoge/real(inum2) ! average
!write(*,*) k,inum2,sla1d(k)
enddo
!----------------------------
! Least-square fitting (1-D)
!----------------------------
B=0; C=0; D=0; E=0;
do k=1,inum
B=B+ここは考えましょう
C=C+ここは考えましょう
D=D+ここは考えましょう
E=E+ここは考えましょう
enddo
slope=ここは考えましょう
offset=ここは考えましょう
C------------------------------
! Build SLA timeseries based on results above
! sla1d2 is the linear trend of SLA
C------------------------------
do k=1,inum
sla1d2(k)=ここは考えましょう
enddo
!-----------------------------------------------
! Deseasonalized but keep the number of data!
!-----------------------------------------------
! 1st step: cal seasonal variation of sla relative to linear trend
! do i=1,12
! icount=0
! monmean(i)=0
! do k=1,inum
! if(ここは考えましょう)then
! monmean(i)=monmean(i)+(ここは考えましょう)
! icount=icount+1
! endif
! enddo
! monmean(i)=monmean(i)/real(icount)
! enddo
! do i=1,12
! write(*,*) "monmean(",i,")=",monmean(i)
! enddo
! 2nd step: erase seasonal trend relative to sla trend
! do i=1,12
! do k=1,inum
! if(ここは考えましょう)then
! sla1df(k)=ここは考えましょう
! endif
! enddo
! enddo
!
!
!----------------------------
! Least-square fitting (1-D) applying to "sla1df"
!----------------------------
! B=0; C=0; D=0; E=0;
! do k=1,inum
! B=B+ここは考えましょう
! C=C+ここは考えましょう
! D=D+ここは考えましょう
! E=E+ここは考えましょう
! enddo
! slope1=ここは考えましょう
! offset1=ここは考えましょう
! write(*,*) "slope1,offset1",slope1*12,offset1
C------------------------------
! Build SLA timeseries based on results above
! sla1df2 is the linear trend of SLA, WITHOUT SEASONAL VARIATION
C------------------------------
! do k=1,inum
! sla1df2(k)=ここは考えましょう
! enddo
!----------------------------
C OUTPUT-->write ascii file
C----------------------------
decim=1./24.!(365/12/2/365)
decim2=decim*2
yr=1993
mo=1
open(80,file='GlobalmeanSLA.txt')
open(81,file='GlobalmeanSLA2.txt')
do k=1,inum
date=yr+decim+decim2*(mo-1)
write(80,100) date,sla1d(k)*100
write(81,100) date,sla1d2(k)*100
mo=mo+1
if(mo.eq.13)then
mo=1
yr=yr+1
endif
enddo
close(80)
close(81)
100 format(2f10.4)
stop
end
地衡流、渦運動エネルギー、渦度の計算
参考資料