Google Graphで棒グラフ

GE-GraphというGoogle Graph上に棒グラフをプロットできるフリーソフトがある。 http://www.sgrillo.net/googleearth/gegraph.htm
これを使ってみた。

結果から。下の図は正味海洋-大気熱フラックス（海から大気）の１月長期平均をGoogle Earth上で棒グラフにしたものである。正値のみ使っている。

このkmzファイルはこちら。

GE-Graphに入力するデータとして
場所名,緯度,経度,値
で構成されるcsvファイルを用意する必要がある。
上の例のcsvファイルはコレ（その１,その２）
GE-Graphが5000地点以上を受け付けないのでファイルを２つに分けている。場所名は特に意味は無いので全て0を入れている。

上のデータをGE-Graphに読み込ませる。設定は以下の図のようにした。

それぞれのデータで"Run"を行ってkmzファイルを作った。Google Earth上で二つのファイルをまとめて、saveしたのが上のkmzファイルである。

最後にcsvファイルを作るpythonプログラムを添付する。
NCEP/NCARのデータをESRLのサイトからPython Interface to GrADS を使って読み込む。
latent/sensible heat flux, longwave/shortwave radiationを全て足して正味の熱フラックスを計算する。
陸・海マスクを読み、海のデータだけを使う。また正（海から大気）のデータだけ残している。データは4000点づつセーブする。

from grads.ganum import GaNum  # for python interface for grads
import numpy as np

#====== DATA read
ga = GaNum(Bin='grads')
ga("sdfopen http://www.esrl.noaa.gov/psd/thredds/dodsC/Datasets/ncep.reanalysis.derived/surface_gauss/lhtfl.mon.ltm.nc") #net latent heat flux
ga("sdfopen http://www.esrl.noaa.gov/psd/thredds/dodsC/Datasets/ncep.reanalysis.derived/surface_gauss/shtfl.mon.ltm.nc") #net sensible heart flux
ga("sdfopen http://www.esrl.noaa.gov/psd/thredds/dodsC/Datasets/ncep.reanalysis.derived/surface_gauss/nswrs.mon.ltm.nc") #net shortwave radiation
ga("sdfopen http://www.esrl.noaa.gov/psd/thredds/dodsC/Datasets/ncep.reanalysis.derived/surface_gauss/nlwrs.mon.ltm.nc") #net longwave radiation
ga("set t 1")  # January
ga("netheat=lhtfl.1+shtfl.2+ nswrs.3+nlwrs.4") # all net heat
netheat=ga.exp("netheat")
ga("close 4")
ga("close 3")
ga("close 2")
ga("close 1")
ga("sdfopen http://www.esrl.noaa.gov/psd/thredds/dodsC/Datasets/ncep.reanalysis/surface_gauss/land.sfc.gauss.nc")  #land sea mask
land=ga.exp("land")
del ga

lon=netheat.grid.lon
lat=netheat.grid.lat

#===== DATA manupilation
data=netheat.data
data[land>0]=-999.  # remove values on land

x,y=np.meshgrid(lon,lat)

vl=0.               # data criteria 

xout=x[data>vl].flatten()      # only  values >vl  are used
xout[xout>180]=xout[xout>180]-360. #   -180.=<lon<=180.  
yout=y[data>vl].flatten()      # only  values >vl  are used
zout=data[data>vl].flatten()   # only  values >vl  are used
points=np.zeros(data[data>vl].size) # position number (anythin is ok)

#==== DATA output
nn=zout.size/4000+1
for n in range(nn):
    st=4000*n
    en=4000*(n+1)
    np.savetxt("data%(n)s.txt" %locals(),np.transpose((points[st:en],yout[st:en],xout[st:en],zout[st:en])),fmt='%4i,%7.3f,%7.3f,%7.3f')