この記事は試行錯誤中の段階かつナンセンスなコードです。自己責任でお願いします。
この記事は古いです。下記にてアップデートしました。
https://dandango.pw/?p=1208
目的
気象庁一元化震源データは未だに固定長データで運用されており、著者は大変ストレスに感じています。
そこで、今回はこの固定長データを、どうにかしてPython3でカンマ区切りにしてみます。
Fortranなら簡単にできるとか言わない。
震源レコードフォーマット
面倒くさいですね。ゴリゴリカンマ区切りに変換してみます。
震源決定方法
A1 震源評価 震源を決定するにあたっての初期条件等
不明の場合は空白
1: 深さフリー
2: 深さ刻み条件(深さを一定の幅で変化させて計算)で最適解を求めた
3: 深さ固定等、人の判断による
4: Depth Phaseを用いた
5: S-Pを用いた
7: 参考（2016年3月まで）
8: 決定不能または不採用
9: 震源固定(最も早い相の情報を検測した観測点の緯度、経度)
M: Matched Filter法を用いた
面倒くさくなる諸悪の根源です。 これをもとに場合分けが必要になります。
震源固定の場合は後々面倒くさいことになります。深さ固定も同じみたいです。
固定長分割
ここら辺を読みました。pandas.read_fwfを使えば行けそう。
その結果がこちら。マグニチュードの分割だけアクロバティックなことをしています。
def kotei(hypo):
# 固定長に分割する
colspecs = [(0,1),(1,5),(5,7),(7,9),(9,11),(11,13),(13,17),(17,21),(21,24),(24,28),(28,32),(32,36),(36,40),(40,44),(44,49),(49,52),(52,53),(53,54), (54,55),(55,57),(57,58),(58,59),(59,60),(60,61),(61,62),(62,63),(63,64),(64,65),(65,68),(68,92),(92,95),(95,96)]
names = ['recode','year','month','day','hour','minute','second','time_SE','latitude_deg','latitude_min','lat_SD','longitude_deg','longitude_min','lon_SD','depth','dep_SD','mag1_1', 'mag1_2','mag1_type','mag2','mag2_type','travel_time','hypo_judge','hypo_type','max_coefficient','scale_victim','tsunami_victim','region_main','region_sub','hypo_locale','obs_number','determination_way']
df = pd.read_fwf(hypo, colspecs=colspecs, names=names)
return dfおまじない
とりあえずFloat型に変換しておきます。
def convert(input):
if input.latitude_deg.dtypes == 'O' :
input['latitude_deg'] = input['latitude_deg'].str.replace(' ', '')
if input.longitude_deg.dtypes == 'O' :
input['longitude_deg'] = input['longitude_deg'].str.replace(' ', '')
if input.depth.dtypes == 'O' :
input['depth'] = input['depth'].str.replace(' ', '')
if input.hypo_judge.dtypes == 'O' :
input.replace({'hypo_judge' : {'M' : 10}}, inplace=True)
input.hypo_judge = input.hypo_judge.astype("float64")
input.latitude_deg = input.latitude_deg.astype("float64")
input.latitude_min = input.latitude_min.astype("float64")
input.longitude_deg = input.longitude_deg.astype("float64")
input.longitude_min = input.longitude_min.astype("float64")
input.depth = input.depth.astype("float64")
return input度分を直す
F4.2 震央の緯度(分)
震源固定の場合は小数点以下空白
場合分けが必要なので面倒くさいですね。
def degree(input) :
# 度分秒を度のみに直す
# 震源固定の場合は小数点以下空白
# 震源固定は震源評価9番。
input.latitude_deg.mask(input.hypo_judge == 9, input.latitude_deg + (input.latitude_min / 60.0), inplace=True)
input.longitude_deg.mask(input.hypo_judge == 9, input.longitude_deg + (input.longitude_min / 60.0), inplace=True)
## それ以外はF4.2なので、
input.latitude_deg.where(input.hypo_judge == 9, input.latitude_deg + (input.latitude_min / 60.0 * 0.01), inplace=True)
input.longitude_deg.where(input.hypo_judge == 9, input.longitude_deg + (input.longitude_min / 60.0 * 0.01), inplace=True)
return inputマグニチュード１を直す
※便宜上、マグニチュード２は無視しています。
マグニチュード F2.1 1 気象庁マグニチュード、
気象庁CMTのモーメントマグニチュードまたはUSGS等が計算したマグニチュード
0未満の場合は以下のように表記する
-0.1: -1 -0.9: -9 -1.0: A0
-1.9: A9 -2.0: B0 -3.0: C0
マグニチュード1が求まらなかった場合は空白(半角スペース×2)
めちゃくちゃ複雑です。上のコードではマグニチュードを桁ごとに分けて処理していますが、自分でも正しいのか怪しいです。
def magnitude(input) :
# マグニチュードを直す
# 0未満の場合は以下のように表記する
# -0.1: -1 -0.9: -9 -1.0: A0
# -1.9: A9 -2.0: B0 -3.0: C0
if input.mag1_1.dtypes == 'O' :
input.replace({'mag1_1' : {'-' : -0}}, inplace=True)
input.replace({'mag1_1' : {'A' : -1}}, inplace=True)
input.replace({'mag1_1' : {'B' : -2}}, inplace=True)
input.replace({'mag1_1' : {'C' : -3}}, inplace=True)
input.mag1_1 = input.mag1_1.astype("float64")
input.mag1_2 = input.mag1_2.astype("float64")
input.mag1_1 = input["mag1_1"].where(input["mag1_1"] >= 0 , input["mag1_1"] - ( input["mag1_2"] / 10.0 ))
input.mag1_1 = input["mag1_1"].where(input["mag1_1"] < 0 , input["mag1_1"] + ( input["mag1_2"] / 10.0 ))
return input発生時刻を直す
14-17 F4.2 秒 オリジンタイムの秒
震源固定の場合は小数点以下空白
これは簡単です。
def minute(input) :
# 秒を直す。
# 震源固定の場合は小数点以下空白
# 震源固定は震源評価9番。
input.second.where(input.hypo_judge == 9, input.second * 0.01, inplace=True )
return input深さを直す
F5.2 深さ(km) 深さフリーの条件で計算した時の震源の深さ(km)
I3, 2X 深さ固定または下記の刻みの条件で計算した時の震源の深さ(km)
10km(1926年～1960年、1967年～1982年)
20km(1961年～1966年)
1km(1983年～)
1982年以前の地震については適宜再調査され、深さフリーまたは1km刻みの震源に置き換えられる
ここは謎です。合っているのか分かりませんが、取り敢えず下記で動いています。
def depth(input) :
# 深さを直す
# 深さフリーならF5.2
# 深さ固定ならI3
# 深さフリーは震源評価1番
input.depth.mask(input.hypo_judge == 1, input.depth * 0.01, inplace=True )
input.depth.mask(input.hypo_judge == 10, input.depth * 0.01, inplace=True )
return inputこれ、
input.depth.where(input.hypo_judge == 3, input.depth * 0.01, inplace=True )
でもいいんですかね？よく分かりません。
整理
あとはお掃除です。
FAR FIELDは整理する際に厄介になるので削除しています。
# FAR FIELDは削除する
df = df[df.hypo_locale != "FAR FIELD"]
# 出力
df2 = df.round(4)
df2 = df2.drop(columns=['latitude_min','longitude_min', 'mag1_2'])
df2.to_csv(str(input) + ".csv", index=None)とすれば、良い感じになると思います。
まとめ
こんな感じで組み合わせればうまく行くと思います。
def main() :
input = args[1]
# 固定長データに分割
df = kotei(input)
# 型を変更する
df = convert(df)
# 度分を直す
df = degree(df)
# マグニチュードを直す
df = magnitude(df)
# 秒を直す
df = minute(df)
# 深さを直す
df = depth(df)
# FAR FIELDは削除する
df = df[df.hypo_locale != "FAR FIELD"]
# 出力
df2 = df.round(4)
df2 = df2.drop(columns=['latitude_min','longitude_min', 'mag1_2'])
#
# df2.to_csv(str(input) + ".csv", index=None, header=None)
df2.to_csv(str(input) + ".csv", index=None)本当に面倒くさいのでカンマ区切りの震源レコードも作成してほしいです。。
USGSは作成しているのに…。
今日はここまで。
コメント