Fehler beim regridding 3-D-Satellitendaten in Python mit Basemap, 2-D arbeiten

Question

Ich habe folgendes Problem:

Ich habe Wolkendecke Datensätze von verschiedenen Satelliten, die ich auf das Gitter eines Klima regrid will Modell, um Vergleiche zwischen der Modellausgabe und den beobachteten Satellitendaten zu ziehen.

Für jetzt verwendete ich die interp-Funktion von Grundkarte, die perfekt für Arrays in der Form funktioniert: 1 x Länge x Breite, aber es funktioniert nicht für Arrays der Form n x Länge x Breite. Was wäre der beste Weg, um diese Art von 3-D-Arrays neu zu regulieren?

from mpl_toolkits.basemap import interp 
import xarray as xr 
def regrid_MAR10km(x_in,y_in,data_in): 
    mar_10km = xr.open_dataset('/media/..../MAR_10km /MARv3.5.2-10km-ERA-2008.nc') 
    lat = mar_10km['LAT'] 
    lon = mar_10km['LON'] 
    result = interp(data_in, x_in, y_in,lon,lat) 
    return result 

Mein Problem ist, dass ich Fehlermeldungen der folgenden Form erhalten, wenn ich versuche, 3-D-Daten zu verwenden, ist das Array in meinem Fall der Form 161 (die die Wolkendecke für jeden Monat ist!) x lon x lat

<xarray.DataArray 'Cloud_Fraction_Mask_Total_Mean' (Begin_Date: 161, lat: 28, lon: 110)> 
[495880 values with dtype=float64] 
Coordinates: 
* Begin_Date (Begin_Date) datetime64[ns] 2002-07-01 2002-08-01 2002-09-01 ... 
* lat   (lat) float32 85.5 84.5 83.5 82.5 81.5 80.5 79.5 78.5 77.5 ... 
* lon   (lon) float32 -94.5 -93.5 -92.5 -91.5 -90.5 -89.5 -88.5 ... 

Und das ist der Fehler, dass es produziert:

IndexError        Traceback (most recent call last) 
<ipython-input-19-5117a62e570e> in <module>() 
----> 1 cloud_regrid =  fc.regrid_MAR10km(longitude,latitude,cloud_data_UD) 

/media/sf_Shared/Black_and_bloom/CODE/functions.py in regrid_MAR10km(x_in, y_in, data_in) 
20   lat = mar_10km['LAT'] 
21   lon = mar_10km['LON'] 
---> 22   result = interp(data_in, x_in, y_in,lon,lat) 
23   return result 

/home/sh16450/miniconda3/envs/snowflakes/lib/python3.5/site-packages/mpl_toolkits/basemap/__init__.py in interp(datain, xin, yin, xout, yout, checkbounds, masked, order) 
4958   dataout = (1.-delx)*(1.-dely)*datain[yi,xi] + \ 
4959     delx*dely*datain[yip1,xip1] + \ 
-> 4960     (1.-delx)*dely*datain[yip1,xi] + \ 
4961     delx*(1.-dely)*datain[yi,xip1] 
4962  elif order == 0: 

IndexError: index 41 is out of bounds for axis 1 with size 28 

Answer 1

Vom Basemap docs die interp Methode befasst sich ausschließlich mit 2-D-Arrays, wobei die erste Dimension der km entspricht imension (Breitengrade) und der 2. bis x (Längengrade). Wenn Sie etwas anderes geben, fürchte ich, dass Sie mit schlechten Ergebnissen enden werden.

Was tatsächlich passiert, ist, dass, da Basemap annimmt, dass Ihre Daten 2D sind, versucht wird, die Interpolation auf Ihren Zeit- und Breiten- (161,28) Achsen so auszuführen, als wären sie Breiten- und Längengrad. Da das Longitude-Array (110,), das Sie an die Methode übergeben haben, jetzt größer ist als die Achse, die Basemap als Longitude (28,) annimmt, haben Sie am Ende eine Formabweichung, die zu diesem IndexError führt.

Ich bin nicht genau sicher, warum die (1,28,110) funktioniert, aber ich wette, dass es einige Array flattening oder squeezing, dass NumPy unter den Abdeckungen, die im Grunde ignoriert diese leere Dimension endet.

Ich empfehle entweder die einzelnen Zeitschritte zu durchlaufen und die Interpolation auf diesem 2D-Datenabschnitt durchzuführen oder die scipy.interpolate Dokumentation zu durchsuchen, um zu sehen, ob sie eine Funktion bieten, die das macht, wonach Sie suchen.

Viel Glück!