MODIS aqua Global Level 3 Mapped Mid-IR 海水面温度をPyNIOで読み込み、PyNGLでプロットした(daily, 4km)。
データはNASA PO.DAACのサイトのFTPサイトより、8月3日のデータを入手し(A2010215.L3m_DAY_SST4_4.bz2)、bunzip2で解凍しておく。
データはhdf4フォーマットである。
以下がスクリプト。
SSTは変数名"l3m_data"でunsigned 16bit integer
クオリティフラッグが変数名"l3m__qual"でunsigned 8bit integer
これらはHDFVIEW で確かめた。
クオリティフラッグが0以外は欠損値とした。
データは、南北が逆になっているので反転してある。
import Nio import Ngl import numpy as np import scikits.timeseries as ts #day year=2010 month=8 day=3 t=ts.Date("D",year=year,month=month,day=day) days='%03d' % t.day_of_year td=str(year)+"/"+str(month)+"/"+str(day) # open file nfile=Nio.open_file("A%(year)s%(days)s.L3m_DAY_SST4_4" %locals(), mode='r',format="h4") nfile.set_option("MaskedArrayMode","MaskedNever") #grid lon=nfile.Westernmost_Longitude+nfile.Longitude_Step[0]*np.arange(nfile.Number_of_Columns) lat=nfile.Southernmost_Latitude+nfile.Latitude_Step[0]*np.arange(nfile.Number_of_Lines) #grid salection xrange=np.logical_and(lon>=130,lon<=140.) yrange=np.logical_and(lat>=30,lat<=36.) temp=nfile.variables["l3m_data"] flag=nfile.variables["l3m_qual"] temp2=((temp[::-1,:])[yrange,:])[:,xrange].astype("uint16") flag2=((flag[::-1,:])[yrange,:])[:,xrange].astype("uint8") sst=temp2*temp.Slope[0]+temp.Intercept[0] sst=np.ma.masked_array(sst,mask=flag2>0) nfile.close() # plot by PyNGL # Open a workstation. # wks_type = "png" wks = Ngl.open_wks(wks_type,"MODIS") resources = Ngl.Resources() resources.tiXAxisString = "~F25~longitude" resources.tiYAxisString = "~F25~latitude" resources.cnFillOn = True # Turn on contour fill. resources.cnLineLabelsOn = False # Turn off line labels. resources.cnInfoLabelOn = False # Turn off info label. resources.mpGridAndLimbOn = False resources.sfXCStartV = float(min(lon[xrange])) # Define where contour plot resources.sfXCEndV = float(max(lon[xrange])) # should lie on the map plot. resources.sfYCStartV = float(min(lat[yrange])) resources.sfYCEndV = float(max(lat[yrange])) resources.mpLimitMode = "LatLon" # Limit the map view. resources.mpMinLonF = float(min(lon[xrange])) resources.mpMaxLonF = float(max(lon[xrange])) resources.mpMinLatF = float(min(lat[yrange])) resources.mpMaxLatF = float(max(lat[yrange])) resources.mpDataBaseVersion = "MediumRes" resources.tiMainString = "MODIS aqua Global Level 3 Mapped Mid-IR SST day=%(td)s" %locals()# Set a title. resources.tiMainFontHeightF=0.015 resources.cnLevelSelectionMode = "ExplicitLevels" # Define own levels. resources.cnLevels = np.arange(20.,30.,0.5) # # draw contours over map. # map = Ngl.contour_map(wks,sst,resources) Ngl.end()MODISにはAquaとTerra、それぞれに海面水温データとしてはMid-IR SSTとthermal IR SSTがある。どのように使いわけたらいいのだろうか?
参照
JAXA MODISページ
http://www.eorc.jaxa.jp/hatoyama/satellite/sendata/modis_j.html
EORC
http://kuroshio.eorc.jaxa.jp/ADEOS/mod_nrt_new/index.html
MODIS wikipedia
http://ja.wikipedia.org/wiki/MODIS
0 件のコメント:
コメントを投稿