MM5修改海溫

MM5修改海溫：

　　在進行地形敏感度數值實驗時，除了將地形移除的無地形實驗之外，也常常會進行將地表分類替換為海洋的Ocean Run，以藉此檢驗地表特性、洋面通量對此實驗之影響。

但是問題來了，地表分類是替換掉了，但是地表溫度呢？

　　以台灣本島這種高聳地形為例，高山地區的地表溫度一定明顯低於洋面溫度，因此在這個假想的Ocean Run中，這塊原本是陸地的「海洋」會有著很不尋常的低溫存在，這時就可以修改海溫分布來讓他看起來「正常」點，然而，溫度分布又該採用哪種分布呢？長期平均？線性內插？甚至是spline內插？這種敏感度實驗的狀況本來就不存在於現實之中，所以選用哪種分布，甚至是不做任何改變都是可以討論的，端看你希望透過這個實驗看出什麼。
以下就說明一下怎麼修改海溫，內容同樣是出自我之前寫給學弟妹的MM5使用手冊。

開始修改可以先看一下這張取自MM5 OnLine Tutorial - SST Note的海溫分布圖：
SST
這裡的「海溫」連陸地上都有值，不過不用擔心，根據線上教學的說明呢，陸上的值沒有任何意義，只是經由MM5自己內插產生而成的。

然而，在表面溫度（SkinTemperature）這個變數裡，則是包含了陸地、海洋的溫度資料，從下圖中可以看出海溫差異不大，但是有表現出陸地的溫度特性來，然而由於海陸的溫差較大，在交界處也就理所當然的會有著相當大的溫度梯度。

在REGRID裡面對這兩個場的處理方式是一樣的，不過稍後在INTERPF程式中就不一樣了。

　　在INTERPF中，如果程式在讀取來自REGIRD的資料時有找到SST的資料，那就會將其存入TSEASFC這個變數裡，然後用他的資料當作水體的表面溫度來製作下邊界場（這個變數是可以隨時間變化的）程式也會依據海陸分布的遮罩自動判斷哪裡要用SST的資料，哪裡要用地表溫度的資料。
　　反之，如果沒有找到SST的資料而只有表面溫度的資料，那就會將此表面溫度所有時間的值做平均，產生一個不隨時間變動的TSEASFC海溫場供模式使用，這個平均的目的就是要移除在平均的海溫場中，陸地因日夜變化而造成的溫度差異。

　　劈哩啪啦說了這麼多，但最關鍵的點就是：「模式會自動依據海陸遮罩，來決定哪個位置要調用海溫資料，哪個位置要調用表面溫度資料」。然而需不需要修改海溫則可以先將MMINPUT的資料抓出來繪圖(SST與TG)，看看模式預處理後的結果如何：
修改前的海溫

上面這兩張圖最明顯的差別除了溫度的不同之外，就是在地表溫度中Domain西北方的中國沿海以及東方的日本琉球群島上有著較高的地表溫度，而Domain中間的低值就是我想要修改的部分。

　　對這些變數有初步的了解之後就可以開始動手啦！海溫資料會存在於初始場（MMINPUT）以及邊界場（BDYOUT）中，修改得由 INTERPF/src/modual_all_io.F 著手，開始之前要先準備好修改好的二進位海溫格點資料，格點數要跟你的domain設定一樣，如果手邊沒有海溫資料的話，也可以先在這個檔案中讓他輸出tseasfc這個變數，再依據這個原本的海溫資料來改成你想要的分布。

　　首先來修改邊界條件，假設手邊已經準備好了一份叫做SSTok.dat的海溫資料檔，請先在modual_all_io.F裡面的add_lbc 這支副程式中（我猜他是 add lower bounday condition 的縮寫，從interpf.f中也可以看到這是用來製作下邊界用的）找到下面這段程式碼：

! Specific things for sst small header.

tseasfc_sh%name = 'SST AVG '
tseasfc_sh%units = 'K '
tseasfc_sh%description = 'Temporal mean of skin temperature '

WRITE ( unit_lowerbc ) tseasfc_sh%num_dims, tseasfc_sh%start_dims, tseasfc_sh%end_dims,
　　　　　tseasfc_sh%xtime , &
　　　　　tseasfc_sh%staggering , tseasfc_sh%ordering ,& tseasfc_sh%current_date , 　　　　　tseasfc_sh%name,& tseasfc_sh%units , tseasfc_sh%description

WRITE ( unit_lowerbc ) tseasfc

這裡的第一個write是寫入檔頭資訊用的，下一個write就是寫入海溫變數tseasfc了，所以只要在第二個write之前把變數tesasfc改掉即可：

integer :: i,j
...
print*,"=== Modifing LOWBDY SST ==="
open(99,file="SSTok.dat",access="direct",form="unformatted",recl=4*ie*je,status="old")
read(99,rec=1) ((tseasfc(i,j),j=1,je),i=1,ie) ! 讀進新的海溫資料，記得要先宣告i,j
close(99)
...
WRITE ( unit_lowerbc ) tseasfc

　　接下來是修改輸入場，由同一個檔案裡的outmodel這支副程式開始，先找到以下程式碼：

WRITE ( immout ) ps0
DO loop = 1 , num_2d
IF ( all_2d(loop)%small_header%ordering(1:2) .EQ. 'YX' ) THEN
WRITE ( immout ) sh_flag
......

這裡的迴圈loop是將all_2d裡面的二為分布場依序寫入，而海溫資料在第16筆（老實說我也忘了當初怎麼測出是第16筆的，可能是將所有description印出來吧）。所以跟處理邊界場一樣，先將準備好的海溫資料讀入之後，在程式寫入這些二維分布場之前將他替換掉就好囉：

integer :: i,j
real,allocatable :: MYtseasfc(:,:)
...
allocate(MYtseasfc(ie,je))
...
WRITE ( immout ) ps0
print*,"=== Modifing LOWBDY SST ==="
open(99,file="SSTok.dat",access="direct",form="unformatted",recl=4*ie*je,status="old")
read(99,rec=1) ((MYtseasfc(i,j),j=1,je),i=1,ie) ! 讀入資料，記得要先宣告i,j以及一個暫存陣列
close(99)
all_2d(16)%array(is:ie,js:je)=tseasfc
DO loop = 1 , num_2d
IF ( all_2d(loop)%small_header%ordering(1:2) .EQ. 'YX' ) THEN
WRITE ( immout ) sh_flag
......

若擔心海溫不是all_2d中第16筆資料的話，也可以在迴圈中加入判斷式比對all_2d(loop)%small_header%name是不是TSEASFC，然後再決定要不要替換掉裡面的資料。

修改結束後請回到INTERPF中重新編譯。以下為修改後的值：

修改後之海溫

　　這裡因為只有修改SST，沒有修改TG，所以在表面溫度的變數中仍然看得到那一塊異常低溫值，不過不用擔心，就如同前面所提到的，模式會自己用海陸遮罩判斷要用哪裡的資料，從模式輸出(MMOUT_DOMAIN1_00)中可以證實這一點。下面這兩張分別是海溫輸出，以及表面溫度的輸出值（請注意陸地區域的溫度）：

而在強度上也有著相當顯著的差異：
最大風速

That's it, have fun!

<-Memo