2010年5月27日木曜日

東北大学 外洋域新世代海面水温



東北大学外洋域新世代海面水温(http://www.ocean.caos.tohoku.ac.jp/~merge/sstbinary/actvalbm.cgi)(分解能1/20度)を読むモジュールである。
ファイル名としてはgzipされたファイルを指定でき、解凍する必要はない。
(ftpのURLも直接指定できるがネットワークの負荷を考えること。)
メインプログラムとして走らせた時は領域一部をmatplotlibでプロットする。

readNGSST.py
import numpy as np
import matplotlib.pyplot as plt

def readNGSST(dir,file):
    """
   read Tohoku univ. New Generation Sea Surface Temp.

   usage:
      lon,lat,sst=readNGSST(dir,file)
     
   input
      dir working dir
      file filename 
           gzipped file is fine.
           URL also works
   output
      lon longitude
      lat latitude
      SST sst data (masked array)          
    """
    # read data
    ds=np.DataSource(dir)
    f=ds.open(file, 'r')
    f.seek(200) # skip header
    sst=np.frombuffer(f.read(),dtype="u1")
    f.close()

    # as masked array
    sst=np.ma.masked_array(sst,mask=np.logical_or(sst==255, sst==0))
    sst=sst.reshape(1000,1000)*0.15-3.

    # grid
    lon=116.+np.arange(1000)*0.05
    lat=13.+np.arange(1000)*0.05

    return lon,lat, sst

############################################
if __name__ == '__main__':

    # example plot

    #read data
    filename="msst2009010112.raw.gz"
    lon,lat,sst=readNGSST(dir=".",file=filename)

    #grid salection
    xrange=np.logical_and(lon>=130,lon<=145.)
    yrange=np.logical_and(lat>=30,lat<=45.)
    # plot
    plt.figure()
    plt.pcolor(lon[xrange],lat[yrange],(sst[yrange,:])[:,xrange]) 
    plt.savefig("ngsst.png")
    plt.show() 

参考
書籍
Matplotlib for Python Developers: Build Remarkable Publication Quality Plots the Easy Way

0 件のコメント: