awk運用三維數組進行插值獲得任意經緯度處的水層沉積層地殼厚度

• 2022 年 6 月 18 日
• 筆記
• shell

awk三維數組與插值

目的：給定經緯度，獲得該點地下的冰層水層沉積層和地殼的厚度
實現：awk一行命令

下載Crust1.0模型

該數據集的詳細介紹見官網.
解壓後有幾個文件:crust1.vp,crust1.vs,crust1.rho,crust1.bnds，分別代表P波速度，S波速度，密度和深度.

數據格式

每個文件共有64800行，9列，每行經緯度如下表所列。crust1.bnds的每行代表某個經緯度所對應的水層，冰層，上中下沉積層和上中下地殼的深度。

awk三維數組存儲

在BEGIN{}里將三維數組初始化，在{}里將數據存儲到三維數組中，在END{}里線性插值

lat=-89.5
lon=-179.5
cat /opt/crust10/crust1.bnds | awk 'BEGIN{
  for (lat_index=0; lat_index<180; lat_index++) {
    for (lon_index=0; lon_index<360; lon_index++) {
      for (i=1; i<=9; i++) {
        value[lat_index,lon_index,i]=0.0
      }
    }
  } 
}{
  lat_index=180-int(NR/360); 
  lon_index=NR%360; 
  for (i=1; i<=9; i++) {
    value[lat_index,lon_index,i]=$i
  }
}END{
  lat_index1=int('"$lat"'+89.5)
  lat_index2=lat_index+1
  lon_index1=int('"$lon"'+179.5)
  lon_index2=lon_index+1
  for (i=1; i<9; i++) {
    lat1='"$lat"'+89.5-lat_index1
    lat2=1.0-lat1
    lon1='"$lon"'+179.5-lon_index1
    lon2=1.0-lon1
    print i,value[lat_index1,lon_index1,i]*lon2*lat2+value[lat_index1,lon_index2,i]*lon1*lat2+value[lat_index2,lon_index1,i]*lon2*lat1+value[lat_index2,lon_index2,i]*lon1*lat1
  }
}' | cat