2017-01-24 16 views
0

グリッドデータ(180度の緯度ボックス×360経度ボックス)があります。海氷データ(0 =氷なし、1 =完全氷)を含むここからダウンロードできます。 https://1drv.ms/u/s!ApWaJaYUUERKjqAP-3K48AuUAFgPXQグリッドデータの面積の計算(度単位の等間隔)

北極氷の広がりを計算しようとしています。つまり、グリッドボックスが15%氷(0.15)を超えている場合、グリッドボックスは氷が覆われているとみなされ、そのグリッドボックスの領域がアレイ。最後に総面積が合計されます。私は必要以上に高い16.08の値を取得しています(データは2014年12月、実際の値は1236です)。

私は、所定の解像度では、例えば、エッジ)が、その多くはない私が使用し

コードは以下の通りです:。!

from netCDF4 import Dataset 
import numpy as np 

Data = Dataset('Ice.nc','r') 

ICE = np.squeeze(np.squeeze(Data.variables['sic'][:])) 
Lat = Data.variables['latitude'][:] 
Lon = Data.variables['longitude'][:] 

Ice_Exten = np.zeros((360,180)) 


for i in range(0,360): 
    for j in range(90,180): #just northern hemisphere 

     if ICE[j,i] > 0.15 and ICE[j,i] <= 1.0: 

      Ice_Exten[i,j] = 12321. * np.cos(np.radians(Lat[j])) 

print np.sum(Ice_Exten)/1e6 

すべてのヘルプははるかに高く評価

+0

グリッドボックスの平均緯度、つまり 'np.cos(np.radians(j-89.5)) 'を使用するとどうなりますか? FWIWを使用すると、緯度を外側ループにして、経度ループの外側の余弦を計算することで、ダブルループをより効率的にすることができます。 –

+0

それを変更するのを忘れて、その行は 'Ice_Exten [i、j] = 12321でなければなりません。* np.cos(np.radians(Lat [j]))'(それは大きな違いはありません) –

答えて

0

。むしろ矛盾は、使用されている陸上/海のマスク、補間によるものです。

0

はなぜ代わりに1にそれをroudingの総氷の面積を追加しません?

from netCDF4 import Dataset 
import numpy as np 

Data = Dataset('Ice.nc','r') 

ICE = np.squeeze(np.squeeze(Data.variables['sic'][:])) 
Lat = Data.variables['latitude'][:] 
Lon = Data.variables['longitude'][:] 

Ice_Exten = np.zeros((360,180)) 


for i in range(0,360): 
    for j in range(90,180): #just northern hemisphere 

     if ICE[j,i] > 0.15 and ICE[j,i] <= 1.0: 

      Ice_Exten[i,j] = 12321. * np.cos(np.radians(j-90.)) * ICE[j,i] 

print np.sum(Ice_Exten)/1e6 

出力:コードとは問題がありませんでした

12.7085786161 
+0

(すなわち、ftp://sidads.colorado.edu/DATASETS/NOAA/G02135/north/monthly/data/N_12_extent_v2.1.csvの最後の列) –

+0

「海の氷の広がりは、海氷の少なくとも一部が存在する海域通常、科学者は氷縁を示すための最小濃度の閾値を定義しており、最も一般的なカットオフは15%であり、科学者は15%のカットオフを使用する。衛星観測と地上観測の間にある。 https://nsidc.org/cryosphere/quickfacts/seaice.html –