2010年3月4日木曜日

串本と浦神の潮位差 (日平均データ)



以前、串本と浦神の潮位差の月平均データを読んだが、今回は日データを読む。
データは気象庁のサイトで手に入れることができる(速報値)。
http://www.data.kishou.go.jp/kaiyou/db/kobe/kuroshio/chouisa/chouisa.php

以下は、データを読むモジュールと時系列プロットするプログラム。

2010/4/6 追記 データのフォーマットと値に修正があり、それに伴いプログラムと図を変更してある。

 

import numpy as np
import scikits.timeseries as ts

def read_kushimoto_uragami_daily():
    day=[]
    value=[]
    ds = np.DataSource(None)
    f=ds.open('http://www.data.kishou.go.jp/kaiyou/db/kobe/kuroshio/chouisa/chou
isa.dat')
    l=0
    for line in f:
         l+=1
         if (l$gt;1): # skip one line
             fields =line.strip().split()
             if len(fields)$gt;=1: day.append(fields[0])
             if len(fields)==1:
                 value.append(-999.)
             elif len(fields)==2:
                 try:
                     a=float(fields[1]) # check if it is number
                     value.append(a)
                 except:
                     value.append(-999.) # if data is not number
    f.close()
    valuea=np.asarray(value,dtype='float32')
    kushimoto_uragami=ts.time_series(valuea,start_date=day[0],freq='d',mask=(val
uea<-990.))    
    return kushimoto_uragami

if __name__ == '__main__':
    import matplotlib.pyplot as plt
    import scikits.timeseries.lib.plotlib as tpl

    kushimoto_uragami=read_kushimoto_uragami_daily()

    fig = tpl.tsfigure()
    fsp = fig.add_tsplot(111)
    fsp.tsplot(kushimoto_uragami, '-')
    plt.xlabel("YEAR")
    plt.ylabel("Kushimoto-Uragami (cm)")
    plt.savefig("kushimoto_uragami_daily.png")
    plt.show()

0 件のコメント: