Mehr

Plotten einer NetCDF-Datei mit lat und lon, die in Variablen enthalten sind

Plotten einer NetCDF-Datei mit lat und lon, die in Variablen enthalten sind


Ich versuche, diese NetCDF-Datei (und ähnliches) als 2D-Konturkarte zu zeichnen. Die Datei enthält eine einzelne Ebene einer Variablen, "pr".
Die Datei hat 413x229 X-Y-Pixelkoordinaten und, gespeichert in den Variablen "xlon" und "xlat", die genauen Koordinaten jedes Gitterpunkts.
Ich kann leicht mit X-Y-Koordinaten plotten:

Dateiname <- "file.nc" Varname <- "pr" Bibliothek(Raster) Variable <- Raster(Dateiname, Varname=Varname) Plot(Variable)

Allerdings kommt der Plot "gequetscht" heraus, und außerdem möchte ich mit Lat-Lon-Koordinaten plotten. Wie kann ich das machen?


Das ist meiner Meinung nach irgendwie ähnlich zu dieser Frage:
https://gis.stackexchange.com/questions/48518/how-to-re-project-ease-equal-area-scalable-earth-grid-with-a-25-km-cylindrica?newreg=0eac93c7cbf04b6386d4e62e64e747ef

Die dort angegebene Lösung ist jedoch für meinen Fall nicht geeignet, weil:

  • Ich möchte vermeiden, die zu verwendenncdf4Paket, würde ich lieber verwendenRasternur weil ich damit das Verfahren auf größere Dateien ausdehnen könnte
  • Ich erhalte mehrere Fehler in derconvertGrid()Funktion, wie zum Beispiel:
    Koordinaten(xp)=~x+y
    Fehler in .subset(x, j) : nur Nullen dürfen mit negativen Indizes gemischt werden

Ich sollte die Projektion umkehren, denke ich. Alles, was ich über die Projektion weiß, ist in der Datei ncdump enthalten:

Projektion = "LAMCON" grid_size_in_meters = 12000. latitude_of_projection_origin = 39. longitude_of_projection_origin = 14. standard_parallel = 35., 51. grid_factor = 0.684241343018562

Die Datei hat keine CRS-Spezifikation. Ich weiß, dass ich in der Lage sein sollte, dies zu tunsp,RasterundrasterVisPakete, aber bisher ist mir das nicht gelungen.

Der Versuch, als data.frame zu zeichnen, funktioniert auch nicht, da das Raster nicht konstant ist.


Ich vermute, dass meine Antwort aus Ihrer anderen StackOverflow-Frage Sie nicht in die richtige Richtung geführt hat?

Hier ist eine detailliertere Antwort, die Sie in die richtige Richtung lenken kann.

Zuerst müssen wir die Projektion Ihrer Daten und die Ausdehnung kennen, um korrekt auf das Long/Lat-Raster projizieren zu können. Leider haben wir das PROJ4 CRS oder den Umfang nicht, daher versuche ich mein Bestes, um diese einzuschätzen.

# Schätzen Sie zunächst das Ausmaß ab. # Wir beginnen mit lat und long und projizieren es dann in die Lambert-konforme Projektionsbibliothek (raster) inputfile <- "pr_averaged_Med_CORDEX_ATM.1980-2008.nc" # Holen Sie sich die lat und lon aus den Daten lat <- raster(inputfile, varname ="xlat") lon <- raster(inputfile, varname="xlon") # Konvertieren Sie in Punkte und passen Sie die Breiten und Längen an plat <- rasterToPoints(lat) plon <- rasterToPoints(lon) lonlat <- cbind(plon[, 3], plat[,3]) # Spezifiziere die lonlat als räumliche Punkte mit Projektion als long/lat lonlat <- SpatialPoints(lonlat, proj4string = CRS("+proj=longlat +datum=WGS84")) # Brauche das rgdal-Paket um es auf die ursprüngliche Koordinatensystembibliothek("rgdal") zu projizieren # Meine beste Vermutung für den proj4-String aus den gegebenen Informationen mycrs <- CRS("+proj=lcc +lat_1=35 +lat_2=51 +lat_0=39 +lon_0 =14 +k=0.684241 +units=m +datum=WGS84 +no_defs") plonlat <- spTransform(lonlat, CRSobj = mycrs) # Schau mal plonlat Extent(plonlat) # Yay! Jetzt können wir die Koordinateninformationen für das Raster richtig einstellen pr <- raster(inputfile, varname="pr") # Korrigieren Sie die Projektion und das Ausmaß look pr plot(pr) # Auf langes Lat-Gitter projizieren r <- projectRaster(pr, crs=CRS("+proj=longlat +datum=WGS84")) # Schau es dir an r plot(r) # Konturen hinzufügen contour(r , add=TRUE) # Länderlinien hinzufügen library("maps") map(add=TRUE, col="blue")

Ich denke, das sieht richtig aus. Die Landesgrenzen in Meeresnähe scheinen mit den Farben der Handlung zu harmonieren. Ich denke, du kannst es von dort nehmen. Sie können schönere Plots erstellen mitrasterVisoderggplot.


Manchmal, wie bei einem Variablenextrakt aus einer HDF-Datei (hdf5), haben Sie nicht den regulären Zellenabstand. Der Lat- und Lon-Auszug aus der Datei zum Beispiel haben unregelmäßige Abstände. Sie werden also Probleme haben, ein korrektes Raster oder Bild zu erstellen. Daher können Sie die Daten interpolieren. Ich verwende dazu die Funktion interp aus dem Paket akima. Auf diese Weise löse ich das Problem, wenn es eine Matriz oder einen Datensatz mit unregelmäßigen Werten von Lat und Lon gibt. Ich denke, dies ist auch für ggplot2 nützlich.


1. Einleitung

NetCDF ist ein weit verbreitetes Format für den Austausch oder die Verteilung von Klimadaten und wurde auch in anderen Bereichen, insbesondere in der Bioinformatik und in anderen Disziplinen, in denen große mehrdimensionale Datenfelder erzeugt werden, übernommen. NetCDF-Dateien sind selbsterklärend, in dem Sinne, dass sie Metadaten enthalten, die beschreiben, was in einer Datei enthalten ist, z , oder Offsets und Skalierungsfaktoren, die möglicherweise zum Komprimieren der Daten verwendet wurden. NetCDF-Dateien sind auch maschinenunabhängig weil zwischen Servern und Computern mit unterschiedlichen Betriebssystemen übertragen werden kann, ohne die Dateien in irgendeiner Weise konvertieren zu müssen. Ursprünglich für die Speicherung und Verteilung von Klimadaten entwickelt, wie sie beispielsweise durch Klimasimulations- oder Reanalysemodelle generiert werden, können das Format und die Protokolle auch für andere gerasterte Datensätze verwendet werden. NetCDF-Bibliotheken werden von Unidata http://www.unidata.ucar.edu/software/netcdf/ entwickelt und gepflegt, und es gibt einfach zu verwendende Anwendungen zum Erstellen einfacher Visualisierungen von NetCDF-Dateien, wie z. B. Panoply, http://www. giss.nasa.gov/tools/panoply/.

Es gibt zwei Versionen von netCDF netCDF3, das weit verbreitet ist, jedoch einige Größen- und Leistungsbeschränkungen aufweist, und netCDF4, das größere Datensätze unterstützt und zusätzliche Funktionen wie die Dateikomprimierung enthält.

R kann netCDF-Dateien lesen und schreiben (und damit analysieren), indem es die von David Pierce bereitgestellten Pakete ncdf und ncdf4 und andere Pakete wie raster , metR1 und RNetCDF verwendet. Die Pakete ncdf4.helpers und easyNCDF stellen einige zusätzliche Werkzeuge bereit.

Das ncdf4-Paket ist sowohl für Windows als auch für Mac OS X (und Linux) verfügbar und unterstützt sowohl das ältere NetCDF3-Format als auch netCDF4. (Siehe die ncdf/ncdf4-Webseite unter http://cirrus.ucsd.edu/


Plotten einer NetCDF-Datei unter Verwendung von Lat und Lon, die in Variablen enthalten sind - Geographische Informationssysteme

Sollten Array-Indizes bei 0 oder 1 beginnen? Mein Kompromiss von 0,5 wurde, wie ich dachte, ohne gebührende Überlegung abgelehnt. --Stan Kelly-Bootle

Einführung

Bibliothekssoftware wie netCDF oder HDF5 bietet Zugriff auf mehrdimensionale Daten über Array-Indizes, aber wir würden oft lieber auf Daten über Koordinaten zugreifen, wie zum Beispiel Punkte auf der Erdoberfläche oder Raum-Zeit-Begrenzungsboxen.

Die Klima- und Prognose-(CF)-Metadatenkonventionen bieten vereinbarte Möglichkeiten zur Darstellung von Koordinatensystemen, gehen jedoch nicht darauf ein, den Zugriff durch Koordinaten praktisch und effizient zu gestalten. Insbesondere kann die Verwendung von Koordinatensystemen, wie sie im Abschnitt "Zweidimensionale Breiten-, Längen- und Koordinatenvariablen" der CF-Konventionen beschrieben sind, aus verschiedenen Gründen keinen direkten Weg zum Auffinden von Feldindizes aus Breiten- und Längenkoordinaten bieten. Beispiele für solche Daten sind:

  • Verwendung von krummlinigen Gittern, die Küstenlinien folgen
  • Satellitenstreifen georeferenziert zu Lat-Lon-Gittern
  • Koordinatensysteme ohne CF-Attribut grid_mapping_name
  • ungegitterte "Stationsdaten, wie Punktbeobachtungen und Sondierungen

Mit einer realen Welt (naja, OK, und idealisiert sphärisch Earth) sehen wir verschiedene Möglichkeiten für den Zugriff auf Daten nach Koordinaten und wie die Verwendung geeigneter Software und Datenstrukturen die Effizienz eines solchen Zugriffs erheblich verbessern kann.

Ein Beispieldatensatz

Der von uns verwendete Beispieldatensatz stammt aus dem HYbrid Coordinate Ocean Model (HYCOM).

Dieser Datensatz ist aus mehreren Gründen ein gutes Beispiel:

  • Es ist typisch für viele Ozean-, Klima- und Vorhersagemodelldatenprodukte mit mehreren Variablen im selben Raster und entspricht den CF-Konventionen zur Darstellung eines Koordinatensystems.
  • Es verfügt über ein nicht triviales Geokoordinatensystem, das durch 2D-Breiten- und Längengrad-Arrays dargestellt wird, die Hilfskoordinaten für andere Variablen sind. Das Gitter scheint eine Teilmenge eines tripolaren Gitters zu sein, aber keiner der Parameter, die zum Ableiten einer Formel erforderlich sind, für die Elemente einer angegebenen Koordinate entsprechen, ist in der Datei enthalten.

Obwohl für das Thema dieses Blogs nicht sehr relevant, nutzt der Beispieldatensatz auch die Komprimierung mit Füllwerten über Land für Ozeandaten. Die netCDF-Classic-Datei ist 102 MB groß, aber die entsprechende netCDF-4-Classic-Modelldatei mit Komprimierung ist nur 30 MB groß. Nichts in diesem Beitrag hat etwas mit dem netCDF-4-Datenmodell zu tun, es gilt alles für das klassische Modell netCDF-3- oder netCDF-4-Daten.

Hier ist eine Ausgabe von "ncdump -h", die auf die Beispieldatei angewendet wird, relevant für das Problem, das wir präsentieren werden:

Das Zeichnen des gesamten Breiten-Längen-Gitters würde zu einem großen Blob von 600.000 Pixeln führen, aber hier ist eine spärlichere Darstellung des Gitters, die jede 25. Linie parallel zu den X- und Y-Achsen zeigt:

Beispiele für Koordinatenabfragen

Die Temperaturvariable hat zusätzliche Längen- und Breitengradkoordinaten, und wir möchten auf Daten zugreifen, die ihren Werten entsprechen, anstatt auf X- und Y-Array-Indizes. Können wir Werte von netCDF-Variablen wie Temperatur oder Salzgehalt effizient bestimmen, beispielsweise an dem Gitterpunkt, der dem 50. nördlichen Breitengrad und -140.

Allgemeiner gesagt, wie können wir effizient

  • den Wert einer Variablen bestimmen, die einem bestimmten Ort und einer bestimmten Zeit am nächsten ist?
  • auf alle Daten in einem Unterraster zugreifen, das innerhalb eines bestimmten Begrenzungsrahmens für Raum-Zeit-Koordinaten definiert ist?

Warum ein iPython-Notebook für das Beispiel verwenden?

In diesem Beispiel wird der Code zum Lesen von netCDF-Daten und zum Mapping vom Koordinatenraum zum Array-Indexraum in einem iPython-Notebook präsentiert. Warum nicht C, Java, Fortran oder eine andere Sprache, die direkt von Unidata unterstützt wird?

Der Punkt ist Klarheit. Dies ist ein ziemlich langer Blog, der wahrscheinlich schon mehrere "TLDR" -Kommentare gesammelt hat. Der begleitende Beispielcode muss kurz, klar und leicht lesbar sein.

Das iPython-Notebook, das aus einer Session zum Lesen von netCDF-Daten im Unidata 2013 TDS-Python Workshop entstand, ist ab dem 2013 Unidata TDS-Python Workshop erhältlich, mit dem Sie Beispiele und Timings auf Ihrem eigenen Computer und sogar mit Ihre eigenen Daten!

Wenn Sie sich für die zweite Methode entscheiden, können Sie die Ergebnisse der Änderung der Beispiele oder des Ersetzens Ihrer eigenen Daten sehen. Um das Notebook verwenden zu können, müssen Python und einige zusätzliche Bibliotheken installiert sein, wie unten auf der README-Seite auf der referenzierten Workshop-Site aufgeführt.

Ein Blick auf die HTML-Version sollte ausreichen, um zu sehen, was vor sich geht.

Vier Ansätze

Im Notebook implementieren und verbessern wir schrittweise mehrere Möglichkeiten, um auf netCDF-Daten basierend auf Koordinaten anstelle von Array-Indizes zuzugreifen:

  • langsam, einfach und manchmal falsch
  • schnell, einfach und manchmal falsch
  • nicht ganz so schnell, aber richtig
  • schnell, korrekt und skalierbar für viele Abfragen an denselben Datenpunkten

Naive, langsame Herangehensweise

Da unser konkretes Problem darin besteht, herauszufinden, welcher von über 600.000 Punkten dem (lat, lon)-Punkt (50, -140) am nächsten liegt, sollten wir zunächst sagen, was mit "close" gemeint ist.

Eine naive Metrik für die Distanz zum Quadrat zwischen (lat,lon) und (lat0, lon0) verwendet einfach die "euklidische Norm", weil sie einfach zu berechnen ist:

Also müssen wir Indizes finden x und ja so dass der Punkt (Breitengrad[ja, x], Längengrad[ja, x]) liegt nahe am gewünschten Punkt (lat0, lon0). Der naive Weg, dies zu tun, besteht darin, zwei verschachtelte Schleifen zu verwenden und den Abstand zum Quadrat für alle über 600.000 Paare von . zu überprüfen ja und x Werte, Punkt für Punkt.

Der Code für diese Version wird als Funktion dargestellt naiv_langsam und einige nachfolgende Zeilen, die die Funktion in unserem expliziten Beispiel aufrufen.

Beachten Sie, dass die Codezellenbeispiele alle zunächst die Beispieldatei öffnen, die zum Finden der gewünschten Werte erforderlichen Daten lesen und die Datei schließen, nachdem genügend Informationen gedruckt wurden, um zu überprüfen, ob wir die richtige Antwort erhalten haben. Das Öffnen der Datei, das Lesen der Daten und das Schließen der Datei wird für jede Funktion wiederholt, sodass die Zellen unabhängig voneinander sind und in beliebiger Reihenfolge ausgeführt oder erneut ausgeführt werden können, ohne dass Fehler auftreten, z Dateien oder schließen Sie dieselbe Datei zweimal.

Schneller mit NumPy

Der erste Ansatz verwendet konventionelle explizite Schleifen, wodurch er Fortran- oder C-Code ähnelt. Es kann um einen Faktor von über 700 beschleunigt werden, indem man einfach NumPy verwendet, eine Python-Bibliothek, die mehrdimensionale Arrays und Array-zu-A-Time-Operationen und -Funktionen unterstützt, die oft explizite Schleifen durch ganze Array-Anweisungen ersetzen. Wir nennen diese Version der Funktion naiv_fast.

Mögen naiv_langsam, leidet es unter Fehlern bei der Verwendung eines "flachen" Maßes für die Nähe.

Tunnelabstand: immer noch schnell, aber korrekter

Ein Problem bei beiden vorherigen Versionen ist, dass sie die Erde als flach behandeln, mit einem Längengrad in der Nähe der Pole, der genauso groß ist wie ein Längengrad am Äquator. Es behandelt auch den Abstand zwischen Punkten an den Kanten, wie (0.0, -179,99) und (0.0, +179,99) als groß, obwohl diese Punkte genau über der Internationalen Datumsgrenze voneinander liegen.

Eine Möglichkeit, diese beiden Probleme zu lösen, besteht darin, eine bessere Metrik zu verwenden, zum Beispiel die Länge eines Tunnels durch die Erde zwischen zwei Punkten als Entfernung, die mit ein wenig Trigonometrie leicht zu berechnen ist.

Dies führt zu einer korrekteren Lösung, die auch dann funktioniert, wenn Längengrade in einem anderen Bereich angegeben werden, z. Wikipedia hat die einfachen trigonometrischen Formeln, die wir verwenden können. Dies gibt uns auch die gleiche Antwort wie der naive Ansatz für unseren Beispielpunkt, vermeidet jedoch falsche Antworten, die wir erhalten würden, wenn wir die Erde in den "naiven"-Ansätzen als flach behandeln würden.

Timing dieses Ansatzes mit der Funktion, die wir aufrufen tunnel_fast, weil wir immer noch schnelle NumPy-Arrays und -Funktionen anstelle von expliziten Schleifen verwenden, ist immer noch etwa 200-mal so schnell wie der loopy-naive Ansatz.

Verwenden von KD-Bäumen für skalierbarere Beschleunigungen

Schließlich verwenden wir eine Datenstruktur, die speziell entwickelt wurde, um schnell den nächsten Punkt in einer großen Menge von Punkten zu einem Abfragepunkt zu finden: den KD-Baum (auch "k-d Baum"). Die Verwendung eines KD-Baums ist ein zweistufiger Prozess.

Zuerst laden Sie den KD-Baum mit dem Satz von Punkten, innerhalb dessen Sie nach einem nächstgelegenen Punkt suchen möchten, normalerweise die Gitterpunkte oder Positionen eines Satzes von ungegitterten Daten. Wie lange dies dauert, hängt von der Anzahl der Punkte sowie deren Dimensionalität ab, entspricht aber der Sortierzeit n Zahlen, nämlich O(N log(N)). Für den Beispieldatensatz mit über 600.000 Punkten im Raster dauert dies auf unserer Laptop-Testplattform weniger als 3 Sekunden, ist aber merklich länger als die Einrichtungszeit für die Suche nach minimaler Tunneldistanz.

Der zweite Schritt stellt einen Abfragepunkt bereit und gibt den nächsten Punkt oder die nächsten Punkte im KD-Baum zum Abfragepunkt zurück, wobei die Definition von "am nächsten" variiert werden kann. Hier verwenden wir 3D-Punkte auf der kugelförmigen Erde mit Einheitsradius. Die Kombination von Setup und Abfrage wird zunächst in a kdtree_fast Funktion.

Obwohl Sie das Setup und die Abfrage nicht zusammen ausführen können, ist die kdtree_fast Abfrage ist deutlich schneller als in der tunnel_fast Suche. Und so kam es dass der kdtree_fast Ansatz skaliert viel besser, wenn wir einen Satz haben, S, von zu suchenden Punkten und viele Abfragepunkte, für die der nächste Punkt (oder die nächsten Punkte) in S.

Zum Beispiel möchten wir möglicherweise Variablenwerte in einem 100 mal 100 (lat,lon) Unterraster der Domäne, und wie wir sehen werden, bietet der KD-Baum eine viel schnellere Lösung als beide tunnel_fast oder naiv_fast Methoden. Wenn sich die Standorte nicht mit Daten aus mehreren Dateien, Zeiten oder Layern ändern, kann die Einrichtungszeit für den KD-Baum für viele Punktabfragen wiederverwendet werden, was erhebliche Effizienzvorteile bietet. Für gängige Geodatensätze auf der Serverseite könnte derselbe KD-Baum für mehrere Datensätze mit denselben Koordinaten wiederverwendet werden.

Die KD-Baum-Datenstruktur kann auch flexibler verwendet werden, um schnelle Abfragen für die M nächsten Punkte zu einem bestimmten Abfragepunkt bereitzustellen, was zum Interpolieren von Werten nützlich ist, anstatt nur den Wert einer Variablen an einem einzelnen nächsten Punkt zu verwenden.

Die Datenstruktur des KD-Baums ist wie ein Schweizer Taschenmesser für Koordinatendaten, in gewisser Weise analog zu regulären Ausdrücken für Strings. Wir haben die Verwendung von KD-Bäumen für 2D-Lat/Lon-Daten und für 3D-Koordinaten von Punkten auf einer Kugel mit einer Tunneldistanzmetrik gezeigt. KD-Bäume können auch mit vertikalen Koordinaten oder Zeit verwendet werden, wenn Sie einen vernünftigen Weg angeben, um "close" zwischen Punkten in der Raumzeit zu definieren.

Implementierungen von KD-Bäumen sind für C, Java, Fortran, C++ und andere Sprachen als Python frei verfügbar. Eine Java-KD-Baum-Implementierung wird in ncWMS verwendet, dem Web Map Service, der in den TDS-Servern enthalten ist, aber an der University of Reading entwickelt wurde. Eine C/C++-Implementierung von KD-Bäumen wird in Fimex verwendet, einer Bibliothek des Norwegischen Meteorologischen Instituts, die die Interpolation und Unterordnung von Geodaten basierend auf dem Unidata Common Data Model unterstützt. Das im iPython-Notebook-Beispiel verwendete Modul scipy.spatial.cKDTree ist eine von Python aufrufbare C-Implementierung, deutlich schneller als die reine Python-Implementierung in scipy.spatial.KDTree.

Der letzte Abschnitt des Notebooks implementiert unsere vier Funktionen als vier objektorientierte Python-Klassen mit dem Einrichtungscode in einem Klassenkonstruktor und einer einzelnen Abfragemethode. Dies macht es einfacher, die Initialisierungs- und Abfrageteile jedes Ansatzes separat zu timen, was zu einfachen Formeln für die benötigte Zeit für N Punkte und einer Vorhersage führt, wann die zusätzliche Einrichtungszeit für kdtree_fast lohnt sich.

Und der Gewinner ist .

Der Umschlag, bitte. Und der Gewinner ist . YMMV!

Die Gründe, warum Ihre Meilenzahl variieren kann, umfassen:

  • Anzahl der Punkte im Suchsatz
  • Dimensionalität der Suchmenge
  • Komplexität der Distanzmetrik
  • Anzahl der Abfragepunkte
  • Anzahl der nächsten Nachbarn benötigt, z.B. zur Interpolation
  • E/A-Faktoren, einschließlich Caching-, Komprimierungs- und Fernzugriffseffekte
  • Berücksichtigung der nicht-sphärischen Erdform

Wir haben unsere Timing-Ergebnisse für die vier Ansätze für das Beispielproblem im Beispieldatensatz am Ende des iPython-Notebooks zusammengefasst. Vielleicht können wir uns in zukünftigen Blogs einige der oben aufgeführten Faktoren ansehen.

Falls Sie das iPython-Notebook nicht öffnen wollten, um die Ergebnisse zu sehen, hier einige Male auf einer schnellen Laptop-Plattform:

Methode Einrichtung (ms) Abfrage (ms) Einrichtung + 10000
Abfragen (Sek.)
Naiv_langsam3.76 7790 77900
Naiv_fast3.82.46 24.6
Tunnel_fast 27.4 5.14 51.4
Kdtree_fast 2520 0.0738 3.3

Feststellen, wann es sich lohnt, einen KD-Baum zu verwenden dieses Beispiel führte zu:

  • kd_tree_fast übertrifft naive_fast über 1050 Abfragen
  • kd_tree_fast übertrifft tunnel_fast über 490 Abfragen

Der Vorteil der Verwendung von KD-Bäumen ist für mehr Such-Sollwerte viel größer, da die KD-Baum-Abfragekomplexität O(log(N)) beträgt, die anderen Algorithmen jedoch O(N) sind, genau wie der Unterschied zwischen der Verwendung einer binären Suche gegen lineare Suche.

Wenn Sie an mehr zu diesem Thema interessiert sind, ist Fast regridding of large, Complex Geospatial Datasets von Jon Blower und A. Clegg ein gutes kurzes Papier über verschiedene Möglichkeiten zum Regriden von Daten mit KD-Bäumen und anderen Methoden.

Sehr schöner Beitrag! Nur ein paar kleine Kommentare:


  1. Das Papier von Blower und Clegg erklärt, warum der "nächste Nachbar"-Punkt nicht immer dasselbe ist wie "die Gitterzelle, die den interessierenden Punkt enthält". Dies kann ein wichtiger Unterschied sein, wenn Sie versuchen, mit verzerrten krummlinigen Gittern genau umzugehen.
  2. Obwohl wir in ncWMS mit einem kd-tree experimentiert haben, fanden wir heraus, dass ein einfacher "Look-Up-Table"-Ansatz schneller war, insbesondere wenn Sie in der Lage sind, die Look-Up-Tabelle zwischenzuspeichern und sie nicht jedes Mal neu generieren zu müssen. Daher verwendet ncWMS standardmäßig die LUT.

Grob gesagt verwendet ncWMS die LUT, um eine Gitterzelle zu finden, die nahe am erforderlichen Lat-Lon-Punkt liegt. Dann sucht es um diesen Punkt herum, bis es die Gitterzelle findet, die diesen Punkt tatsächlich enthält (mit der Java-Methode polygon.contains(point)). Die Methode contains() kann teuer sein, daher ist es im Allgemeinen schneller, den nächsten Nachbarn zu finden, wenn dies wirklich das ist, was Sie wollen.

Ich liebe übrigens die Verwendung von iPython-Notebook!

Geschrieben von Jon Blower am 16. September 2013 um 03:28 MDT #

Danke für die Kommentare, Jon.

Da ich Ihren LUT-Ansatz gelesen habe, möchte ich ihn in demselben Python-Framework implementieren, das ich verwendet habe, und vielleicht einen weiteren Blogeintrag über die Ergebnisse für dieses Beispiel schreiben.

Geschrieben von 67.190.4.13 am 16. September 2013 um 06:21 MDT #

Toll, es wäre sehr interessant, eine Python-Implementierung zu haben. Übrigens hat ncWMS früher einen eher Brute-Force-Algorithmus verwendet, um die LUT zu generieren, was ziemlich langsam war. Es verwendet jetzt eine viel schnellere Methode, bei der die Gitterzellen mit der Java 2D-Grafik-API auf die LUT "gemalt" werden, wodurch die LUT effektiv als Bild behandelt wird. (In der Terminologie des Blower- und Clegg-Papiers wird die LUT mit einem "Source-Push"-Algorithmus generiert.) Sie können den Code unter http://www.resc.rdg.ac.uk/trac/ncWMS/browser einsehen /trunk/src/java/uk/ac/rdg/resc/edal/cdm/LookUpTable.java.</p>

Ich weiß nicht, ob Python äquivalente APIs hat - wenn nicht, könnte der Code ziemlich schwer zu portieren sein!

Geschrieben von Jon Blower am 17. September 2013 um 01:20 MDT #

Ich erinnere mich noch an mein erstes Projekt mit CI und Google Map

und Bibliothekskarten von Google bringen wirklich große Hilfe und einfache Konfiguration

Geschrieben von Qwords Webhosting am 30. November 2015 um 02:43 MST #

Ich stimme zu, dass ich es in demselben Python-Framework implementieren möchte, das ich verwendet habe, und vielleicht einen weiteren Blogeintrag über die Ergebnisse für dieses Beispiel schreiben.

Geschrieben von jasa bikin web am 28. Mai 2017 um 04:08 MDT #

Sie können den Netznummernrechner verwenden. Es ist eine API in Netcdf-Extractor Version 2.0. GNC ist kostenlos und Sie können es nutzen unter:
https://www.agrimetsoft.com/Netcdf-Extractor.aspx</p></p>

Geschrieben von Hasan Sheidaee am 28. April 2018 um 17:29 MDT #

Sie können den Netznummernrechner verwenden. Es ist eine API in Netcdf-Extractor Version 2.0. GNC ist kostenlos und Sie können es nutzen unter:
https://agrimetsoft.com/Netcdf-Extractor.aspx</p></p>

Geschrieben von Hasan Sheidaee am 11. Juli 2018 um 09:34 MDT #

Geschrieben von Hasan shedaee am 17. Juli 2018 um 10:09 MDT #

Hallo! Danke für die Hilfe. Ich kann nicht auf das Notebook zugreifen. Es wurde entfernt. Gibt es eine Möglichkeit auf den Code zuzugreifen. Vielen Dank :)

Geschrieben von Patricio Fernandez am 19. Oktober 2018 um 05:44 MDT #

"Open NC File" kann Daten öffnen und lesen, indem die Liste der Koordinaten von Stationen verwendet wird. Der nc-Benutzer kann die Stationsliste manuell oder per Datei eingeben. "Öffne NC-Datei" kann die Liste der Stationen (bestehend aus Koordinaten) in einer Textdatei oder Excel-Datei oder CSV-Datei lesen. Der nc-Benutzer kann alle nc-Dateien auswählen und in "Öffne NC-Datei" zusammenführen.

Geschrieben von Sohrab Kolsoumi am 19. Oktober 2018 um 08:01 MDT #

Immer noch verwirrend für mich über diesen Programmieralgorithmus. Aber danke fürs teilen

Geschrieben von Jasa-Website Semarang am 14. Dezember 2018 um 19:49 MST #

Hallo! Danke für die Hilfe. Ich kann nicht auf das Notebook zugreifen. Es wurde entfernt. Gibt es eine Möglichkeit auf den Code zuzugreifen. Vielen Dank :)

Geschrieben von Harga Jasa Pembuatan-Website am 16. März 2019 um 04:29 MDT #


Benchmark-Plotfunktionen¶

gcpy.benchmark enthält mehrere Funktionen zum Plotten von GEOS-Chem-Ausgaben in Formaten, die vom GEOS-Chem-Lenkungsausschuss angefordert wurden. Die primäre Verwendung dieser Funktionen besteht darin, Diagramme der meisten GEOS-Chem-Ausgabevariablen zu erstellen, die in bestimmte Kategorien unterteilt sind, z. Spezieskategorien wie Aerosole oder Brom für die SpeciesConc-Diagnostik. In jeder Kategorie erstellen diese Funktionen einstufige PDFs für die Oberfläche und 500hPa sowie PDFs für den zonalen Mittelwert für die gesamte Atmosphäre und nur die Stratosphäre (definiert als 1-100hPa). Für make_benchmark_emis_plots() werden nur Single-Level-Plots an der Oberfläche erzeugt. Alle diese Plotfunktionen enthalten Lesezeichen in den generierten PDFs, die auf die Seiten verweisen, die jede geplottete Menge enthalten. Somit dienen diese Funktionen als Werkzeuge, um schnell umfassende Plots zu erstellen, die zwei GEOS-Chem-Läufe vergleichen. Mit diesen Funktionen werden die öffentlich verfügbaren Plots für 1-Monats- und 1-Jahres-Benchmarks neuer Versionen von GEOS-Chem erstellt.

Viele dieser Funktionen verwenden vordefinierte Variablenlisten (über YAML-Dateien, die in GCPy enthalten sind). Wenn ein Datensatz eine Variable enthält, der andere Datensatz jedoch nicht, werden die Daten für diese Variable im letztgenannten Datensatz als NaN betrachtet und als solche dargestellt.

Gemeinsame Struktur¶

Jede der gcpy.benchmark.make_benchmark_*_plots()-Funktionen erfordert 4 Argumente, um die ref- und dev-Datensätze anzugeben.

Argumente:¶

Pfadname für den Datensatz „Ref“ (auch bekannt als „Referenz“).

Eine Zeichenfolge zur Beschreibung von ref (z. B. Versionsnummer)

Pfadname für den Datensatz „Dev“ (auch bekannt als „Development“).

Dieser Datensatz wird mit dem Datensatz „Referenz“ verglichen.

Eine Zeichenfolge zur Beschreibung von dev (z. B. Versionsnummer)

Beachten Sie, dass die ref- und dev-Argumente in make_benchmark_*_plots() die Pfade zu NetCDF-Dateien sind, anstatt xarray-Datensätze wie in Compare_single_level() und Compare_zonal_mean() . Die Funktionen make_benchmark_*_plots() öffnen diese Dateien intern als Xarray-Datasets und übergeben diese Datasets an Compare_single_level() und Compare_zonal_mean() .

Die Benchmark-Plotfunktionen haben mehrere Schlüsselwortargumente gemeinsam. Schlüsselwortargumente, die nicht den gleichen Zweck für alle Benchmark-Plotfunktionen haben, haben HINWEIS: in der Beschreibung.

Gemeinsame Schlüsselwortargumente:¶

Eine Zeichenfolge, die den Zielordner angibt, in den eine PDF-Datei mit Plots geschrieben wird.

Eine Zeichenfolge, die das Unterverzeichnis von dst angibt, in das PDF-Dateien mit Plots geschrieben werden. In der Praxis wird subdst nur für die 1-Jahres-Benchmarkausgabe benötigt und bezeichnet eine Datumszeichenfolge (z. B. „Jan2016“), die dem Monat entspricht, der geplottet wird. HINWEIS: Nicht verfügbar in wetdep_plots

Setzen Sie dieses Flag auf True, um zuvor erstellte Dateien im Zielordner (angegeben durch das Argument dst) zu überschreiben.

Setzen Sie dieses Flag auf True, um zusätzliche Informationsausgaben zu drucken.

Setzen Sie dieses Flag auf True, um das Plotten von Daten (die oberen beiden Felder jedes Plots, nicht Diffs) auf einer Log-Farbskala zu aktivieren.

sigdiff_files : Liste von str ¶

Dateinamen, die die Liste der Größen enthalten, die signifikante Unterschiede zwischen den Datensätzen aufweisen. Es werden drei Dateien verwendet: eine für die Oberfläche, eine für 500 hPa und eine für den zonalen Mittelwert. Diese Listen werden zum Ausfüllen der Benchmark-Genehmigungsformulare benötigt. HINWEIS: Nicht verfügbar in wetdep_plots

Verzeichnis, das die Datei Spezies_database.yml enthält. Diese Datei wird für Einheitenumrechnungen in ug/m3 verwendet. Die Ausführungsverzeichnisse von GEOS-Chem enthalten eine Kopie dieser Datei, die möglicherweise aktueller ist als die in GCPy enthaltene Version.

Standardwert: Pfad des GCPy-Code-Repositorys

weightsdir : str ¶ Verzeichnis, in dem xESMF-Regridder-NetCDF-Dateien abgelegt (und möglicherweise wiederverwendet) werden. ¶

Definiert die Anzahl gleichzeitiger Arbeiter für paralleles Plotten. Auf 1 setzen, um das parallele Plotten zu deaktivieren. Der Wert -1 ermöglicht der Anwendung die Entscheidung. HINWEIS: In make_benchmark_conc_plots() erfolgt die Parallelisierung auf Artenkategorieebene. In allen anderen Funktionen erfolgt die Parallelisierung innerhalb von Aufrufen von Compare_single_level() und Compare_zonal_mean().

Standardwert: -1 in make_benchmark_conc_plots, 1 in allen anderen

Make_benchmark_aod_plots¶

Funktionsspezifische Schlüsselwort-Argumente:¶

Liste der zu zeichnenden AOD-Variablen. Wenn nicht übergeben, werden alle AOD-Variablen, die Dev und Ref gemeinsam sind, geplottet. Verwenden Sie das Argument varlist, um die Anzahl der Variablen einzuschränken, die beim Debuggen in die PDF-Datei gezeichnet werden.

Diese Funktion erstellt optische Tiefendiagramme der Säulen unter Verwendung des Aerosol-Diagnoseausgangs.

Make_benchmark_conc_plots¶

Funktionsspezifische Schlüsselwort-Argumente:¶

Name der Sammlung, die zum Plotten verwendet werden soll.

Standardwert: „SpeciesConc“

Eine Zeichenfolge, die den Typ der zu zeichnenden Benchmark-Ausgabe angibt, entweder FullChemBenchmark oder TransportTracersBenchmark.

Standardwert: „FullChemBenchmark“

Setzen Sie dieses Flag auf False, um Plots an eine Datei zu senden, anstatt die Datei pro Kategorie zu trennen.

restriktive_katzen : Liste von str ¶

Liste der Benchmark-Kategorien in Benchmark_categories.yml, für die Plots erstellt werden sollen. Wenn leer, werden Diagramme für alle Kategorien erstellt.

Liste der zu erstellenden Plottypen.

Standardwert: [‘sfc’, ‘500hpa’, ‘zonalmean’]

Setzen Sie dieses Flag auf true, um die Normalisierung der Daten nach Fläche zu aktivieren (d. h. kg s-1 –> kg s-1 m-2).

Liste der Benchmark-Kategorien, die in ug/m3 umgerechnet werden sollen

Standardwert: [„Aerosole“, „Secondary_Organic_Aerosols“]

Bereiche: dict of xarray DataArray: ¶

Rasterfeldflächen in m2 auf Ref- und Dev-Rastern.

Pfadname für Referenzmeteorologie

Pfadname für Entwicklungsmeteorologie

Pfadname für einen zweiten „Ref“ (auch „Referenz“) Datensatz für die Diff-of-Diffs-Darstellung. Dieses Dataset sollte denselben Modelltyp und dasselbe Raster wie Ref. aufweisen.

Pfadname für einen zweiten „Ref“ (auch „Referenz“) Datensatz für die Diff-of-Diffs-Darstellung. Dieses Dataset sollte denselben Modelltyp und dasselbe Raster wie Ref. aufweisen.

Diese Funktion erstellt standardmäßig die SpeciesConc-Diagnoseausgabe mit der SpeciesConc-Diagnoseausgabe. Diese Funktion ist die einzige Benchmark-Plotfunktion, die das Plotten von Diff-of-Diffs unterstützt, bei der 4 Datensätze übergeben werden und die Unterschiede zwischen zwei Gruppen von Ref-Datensätzen vs . modellversionsübergreifende Änderungen in GEOS-Chem Classic). Dies ist auch die einzige Benchmark-Plotting-Funktion, die Plots basierend auf der Kategorie (gekennzeichnet durch das plot_by_spc_cat-Flag) an separate Ordner sendet. Die vollständige Liste der Artenkategorien ist in benchmark_categories.yml (in GCPy enthalten) wie folgt aufgeführt:

Make_benchmark_emis_plots¶

Funktionsspezifische Schlüsselwort-Argumente:¶

nach den Benchmark-Spezieskategorien (z. B. Oxidationsmittel, Aerosole, Stickstoff usw.) Diese Kategorien sind in der YAML-Datei Benchmark_species.yml angegeben.

Setzen Sie dieses Flag auf True, um Plots nach HEMCO-Emissionskategorien (z. B. Anthro, Aircraft, Bioburn usw.) in PDF-Dateien aufzuteilen.

Setzen Sie dieses Flag auf True, um die vertikale Ebenenreihenfolge im Datensatz "Ref" umzukehren (falls "Ref" von der Oberseite der Atmosphäre statt von der Oberfläche aus beginnt).

Setzen Sie dieses Flag auf True, um die vertikale Ebenenreihenfolge im Datensatz „Dev“ umzukehren (falls „Dev“ von der Oberseite der Atmosphäre statt von der Oberfläche aus beginnt).

Diese Funktion erzeugt Diagramme der Gesamtemissionen unter Verwendung der Ausgabe von HEMCO_diagnostics (für GEOS-Chem Classic) und/oder GCHP.Emissions-Ausgabedateien.

Make_benchmark_jvalue_plots¶

Funktionsspezifische Schlüsselwort-Argumente:¶

Liste der darzustellenden J-Wert-Variablen. Wenn nicht übergeben, werden alle J-Wert-Variablen, die sowohl dev als auch ref gemeinsam haben, geplottet. Das Argument varlist kann eine nützliche Möglichkeit sein, die Anzahl der Variablen einzuschränken, die beim Debuggen in die PDF-Datei gezeichnet werden.

Setzen Sie dieses Flag, um lokale Mittags-J-Werte darzustellen. Dadurch werden alle J-Wert-Variablen durch den JNoonFrac-Zähler dividiert, der den Bruchteil der Zeit darstellt, in der es an jedem Standort lokaler Mittag war.

Liste der zu erstellenden Plottypen.

Standardwert: [‘sfc’, ‘500hpa’, ‘zonalmean’]

Setzen Sie dieses Flag auf True, um die vertikale Ebenenreihenfolge im Datensatz "Ref" umzukehren (falls "Ref" von der Oberseite der Atmosphäre statt von der Oberfläche aus beginnt).

Setzen Sie dieses Flag auf True, um die vertikale Ebenenreihenfolge im Datensatz „Dev“ umzukehren (falls „Dev“ von der Oberseite der Atmosphäre statt von der Oberfläche aus beginnt).

Diese Funktion erzeugt Diagramme von J-Werten unter Verwendung der JValues ​​GEOS-Chem-Ausgabedateien.

Make_benchmark_wetdep_plots¶

Function-specific keyword args:¶

A string with date information to be included in both the plot pdf filename and as a destination folder subdirectory for writing plots

A string denoting the type of benchmark output to plot, either FullChemBenchmark or TransportTracersBenchmark.


Interpolation Routines¶

Interpolating to a Horizontal Level¶

The wrf.interplevel() function is used to interpolate a 3D field to a specific horizontal level, usually pressure or height.

Vertical Cross Sections¶

The wrf.vertcross() function is used to create vertical cross sections. To define a cross section, a start point and an end point needs to be specified. Alternatively, a pivot point and an angle may be used. The start point, end point, and pivot point are specified using a wrf.CoordPair object, and coordinates can either be in grid (x,y) coordinates or (latitude,longitude) coordinates. When using (latitude,longitude) coordinates, a NetCDF file object or a wrf.WrfProj object must be provided.

The vertical levels can also be specified using the Ebenen Parameter. If not specified, then approximately 100 levels will be chosen in 1% increments.

Example Using Start Point and End Point¶

Example Using Pivot Point and Angle¶

Example Using Lat/Lon Coordinates¶

Example Using Specified Vertical Levels¶

Interpolating Two-Dimensional Fields to a Line¶

Two-dimensional fields can be interpolated along a line, in a manner similar to the vertical cross section (see Vertical Cross Sections ), using the wrf.interpline() function. To define the line to interpolate along, a start point and an end point needs to be specified. Alternatively, a pivot point and an angle may be used. The start point, end point, and pivot point are specified using a wrf.CoordPair object, and coordinates can either be in grid (x,y) coordinates or (latitude,longitude) coordinates. When using (latitude,longitude) coordinates, a NetCDF file object or a wrf.WrfProj object must also be provided.

Example Using Start Point and End Point¶

Example Using Pivot Point and Angle¶

Example Using Lat/Lon Coordinates¶

Interpolating a 3D Field to a Surface Type¶

The wrf.vinterp() is used to interpolate a field to a type of surface. The available surfaces are pressure, geopotential height, theta, and theta-e. The surface levels to interpolate also need to be specified.


Preparation¶

Basic knowledge of Python is assumed. For Python beginners, I recommend Python Data Science Handbook von Jake VanderPlas, which is free available online. Skim through Chapter 1, 2 and 4, then you will be good! (Chapter 3 and 5 are important for general data science but not particularly necessary for working with gridded model data.)

It is crucial to pick up the correct tutorial when learning Python, because Python has a wide range of applications other than just science. For scientific research, you should only read data science tutorials. Other famous books such as Fluent Python and The Hitchhiker’s Guide to Python are great for general software development, but not quite useful for scientific research.

Once you manage to use Python in Jupyter Notebook, follow this GCPy page to set up Python environment for Geosciences. Now you should have everything necessary for working with GEOS-Chem data and most of earth science data in general.


5. Coordinate Systems

A data variable’s dimensions are used to locate data values in time and space or as a function of other independent variables. This is accomplished by associating these dimensions with the relevant set of latitude, longitude, vertical, time and any non-spatiotemporal coordinates. This section presents two methods for making that association: the use of coordinate variables, and the use of auxiliary coordinate variables.

Any of a variable’s dimensions that is an independently varying latitude, longitude, vertical, or time dimension (see Section 1.2, "Terminology") and that has a size greater than one must have a corresponding coordinate variable, i.e., a one-dimensional variable with the same name as the dimension (see examples in Chapter 4, Coordinate Types). This is the only method of associating dimensions with coordinates that is supported by [COARDS].

Any longitude, latitude, vertical or time coordinate which depends on more than one spatiotemporal dimension must be identified by the coordinates attribute of the data variable. Der Wert der coordinates Attribut ist a blank separated list of the names of auxiliary coordinate variables. There is no restriction on the order in which the auxiliary coordinate variables appear in the coordinates attribute string. The dimensions of an auxiliary coordinate variable must be a subset of the dimensions of the variable with which the coordinate is associated, with two exceptions. First, string-valued coordinates (Section 6.1, "Labels") will have a dimension for maximum string length if the coordinate variable has a type of char rather than a type of Schnur . Second, in the ragged array representations of data (Chapter 9, Discrete Sampling Geometries), special methods are needed to connect the data and coordinates

We recommend that the name of a multidimensional coordinate variable should not match the name of any of its dimensions because that precludes supplying a coordinate variable for the dimension. This practice also avoids potential bugs in applications that determine coordinate variables by only checking for a name match between a dimension and a variable and not checking that the variable is one dimensional.

If the longitude, latitude, vertical or time coordinate is multi-valued, varies in only one dimension, and varies independently of other spatiotemporal coordinates, it is not permitted to store it as an auxiliary coordinate variable. This is both to enhance conformance to COARDS and to facilitate the use of generic applications that recognize the NUG convention for coordinate variables. An application that is trying to find the latitude coordinate of a variable should always look first to see if any of the variable’s dimensions correspond to a latitude coordinate variable. If the latitude coordinate is not found this way, then the auxiliary coordinate variables listed by the coordinates attribute should be checked. Note that it is permissible, but optional, to list coordinate variables as well as auxiliary coordinate variables in the coordinates Attribut. If the longitude, latitude, vertical or time coordinate is single-valued, it may be stored either as a coordinate variable with a dimension of size one, or as a scalar coordinate variable (Section 5.7, "Scalar Coordinate Variables").

Wenn ein Achse attribute is attached to an auxiliary coordinate variable, it can be used by applications in the same way the Achse attribute attached to a coordinate variable is used. However, it is not permissible for a data variable to have both a coordinate variable and an auxiliary coordinate variable, or more than one of either type of variable, having an Achse attribute with any given value e.g. there must be no more than one Achse attribute for x for any data variable. Note that if the Achse attribute is not specified for an auxiliary coordinate variable, it may still be possible to determine if it is a spatiotemporal dimension from its own units or standard_name, or from the units and standard_name of the coordinate variable corresponding to its dimensions (see Chapter 4, Coordinate Types). For instance, auxiliary coordinate variables which lie on the horizontal surface can be identified as such by their dimensions being horizontal. Horizontal dimensions are those whose coordinate variables have an Achse attribute of x oder Ja , oder ein units attribute indicating latitude and longitude.

To geo-reference data horizontally with respect to the Earth, a grid mapping variable may be provided by the data variable, using the grid_mapping Attribut. If the coordinate variables for a horizontal grid are not longitude and latitude, then a grid_mapping variable provides the information required to derive longitude and latitude values for each grid location. Wenn nein grid mapping variable is referenced by a data variable, then longitude and latitude coordinate values shall be supplied in addition to the required coordinates. For example, the Cartesian coordinates of a map projection may be supplied as coordinate variables and, in addition, two-dimensional latitude and longitude variables may be supplied via the coordinates attribute on a data variable. The use of the Achse attribute with values x und Ja is recommended for the coordinate variables (see Chapter 4, Coordinate Types).

It is sometimes not practical to specify the latitude-longitude location of data which is representative of geographic regions with complex boundaries. For this purpose, provision is made in Section 6.1.1, "Geographic Regions" for indicating the region by a standardized name.

5.1. Independent Latitude, Longitude, Vertical, and Time Axes

When each of a variable’s spatiotemporal dimensions is a latitude, longitude, vertical, or time dimension, then each axis is identified by a coordinate variable.

xwind(n,k,j,i) is associated with the coordinate values lon(i) , lat(j) , pres(k) , und time(n) .

5.2. Two-Dimensional Latitude, Longitude, Coordinate Variables

The latitude and longitude coordinates of a horizontal grid that was not defined as a Cartesian product of latitude and longitude axes, can sometimes be represented using two-dimensional coordinate variables. These variables are identified as coordinates by use of the coordinates Attribut.

T(k,j,i) is associated with the coordinate values lon(j,i) , lat(j,i) , und lev(k) . The vertical coordinate is represented by the coordinate variable lev(lev) and the latitude and longitude coordinates are represented by the auxiliary coordinate variables lat(yc,xc) und lon(yc,xc) which are identified by the coordinates Attribut.

Note that coordinate variables are also defined for the xc und yc Maße. This faciliates processing of this data by generic applications that don’t recognize the multidimensional latitude and longitude coordinates.

5.3. Reduced Horizontal Grid

A "reduced" longitude-latitude grid is one in which the points are arranged along constant latitude lines with the number of points on a latitude line decreasing toward the poles. Storing this type of gridded data in two-dimensional arrays wastes space, and results in the presence of missing values in the 2D coordinate variables. We recommend that this type of gridded data be stored using the compression scheme described in Section 8.2, "Compression by Gathering". Compression by gathering preserves structure by storing a set of indices that allows an application to easily scatter the compressed data back to two-dimensional arrays. The compressed latitude and longitude auxiliary coordinate variables are identified by the coordinates Attribut.

PS(n) is associated with the coordinate values lon(n) , lat(n) . Compressed grid index (n) would be assigned to 2D index (j,i) (C index conventions) where

Notice that even if an application does not recognize the Kompresse attribute, the grids stored in this format can still be handled, by an application that recognizes the coordinates Attribut.

5.4. Timeseries of Station Data

This section has been superseded by the treatment of time series as a type of discrete sampling geometry in Chapter 9.

5.5. Trajectories

This section has been superseded by the treatment of time series as a type of discrete sampling geometry in Chapter 9.

5.6. Horizontal Coordinate Reference Systems, Grid Mappings, and Projections

EIN grid mapping variable may be referenced by a data variable in order to explicitly declare the coordinate reference system (CRS) used for the horizontal spatial coordinate values. For example, if the horizontal spatial coordinates are latitude and longitude, the grid mapping variable can be used to declare the figure of the earth (WGS84 ellipsoid, sphere, etc.) they are based on. If the horizontal spatial coordinates are easting and northing in a map projection, the grid mapping variable declares the map projection CRS used and provides the information needed to calculate latitude and longitude from easting and northing.

When the horizontal spatial coordinate variables are not longitude and latitude, it is required that further information is provided to geo-locate the horizontal position. EIN grid mapping variable provides this information.

Wenn nein grid mapping variable is provided and the coordinate variables for a horizontal grid are not longitude and latitude, then it is required that the latitude and longitude coordinates are supplied via the coordinates Attribut. Such coordinates may be provided in addition to the provision of a grid mapping variable, but that is not required.

A grid mapping variable provides the description of the mapping via a collection of attached attributes. It is of arbitrary type since it contains no data. Its purpose is to act as a container for the attributes that define the mapping. The one attribute that all grid mapping variables must have is grid_mapping_name, which takes a string value that contains the mapping’s name. The other attributes that define a specific mapping depend on the value of grid_mapping_name. The valid values of grid_mapping_name along with the attributes that provide specific map parameter values are described in Appendix F, Grid Mappings

The grid mapping variables are associated with the data and coordinate variables by the grid_mapping Attribut. This attribute is attached to data variables so that variables with different mappings may be present in a single file. The attribute takes a string value with two possible formats. In the first format, it is a single word, which names a grid mapping variable. In the second format, it is a blank-separated list of words "<gridMappingVariable>: <coordinatesVariable> [<coordinatesVariable> …​] [<gridMappingVariable>: <coordinatesVariable>…​]" , which identifies one or more grid mapping variables, and with each grid mapping associates one or more coordinatesVariables, i.e. coordinate variables or auxiliary coordinate variables.

Where an extended "<gridMappingVariable>: <coordinatesVariable> [<coordinatesVariable>]" entity is defined, then the order of the <coordinatesVariable> references within the definition provides an explicit order for these coordinate value variables, which is used if they are to be combined into individual coordinate tuples.

This order is only significant if crs_wkt is also specified within the referenced grid mapping variable. Explicit 'axis order' is important when the grid_mapping_variable contains an attribute crs_wkt as it is mandated by the OGC CRS-WKT standard that coordinate tuples with correct axis order are provided as part of the reference to a Coordinate Reference System.

Using the simple form, where the grid_mapping attribute is only the name of a grid mapping variable, 2D latitude and longitude coordinates for a projected coordinate reference system use the same geographic coordinate reference system (ellipsoid and prime meridian) as the projection is projected from.

The grid_mapping variable may identify datums (such as the reference ellipsoid, the geoid or the prime meridian) for horizontal or vertical coordinates. Therefore a grid mapping variable may be needed when the coordinate variables for a horizontal grid are longitude and latitude. The grid_mapping_name of latitude_longitude should be used in this case.

The expanded form of the grid_mapping attribute is required if one wants to store coordinate information for more than one coordinate reference system. In this case each coordinate or auxiliary coordinate is defined explicitly with respect to no more than one grid_mapping Variable. This syntax may be used to explicitly link coordinates and grid mapping variables where only one coordinate reference system is used. In this case, all coordinates and auxiliary coordinates of the data variable not named in the grid_mapping attribute are unrelated to any grid mapping variable. All coordinate names listed in the grid_mapping attribute must be coordinate variables or auxiliary coordinates of the data variable.

In order to make use of a grid mapping to directly calculate latitude and longitude values it is necessary to associate the coordinate variables with the independent variables of the mapping. This is done by assigning a Standardname to the coordinate variable. The appropriate values of the Standardname depend on the grid mapping and are given in Appendix F, Grid Mappings.

A CF compliant application can determine that rlon and rlat are longitude and latitude values in the rotated grid by recognizing the standard names grid_longitude und grid_latitude . Note that the units of the rotated longitude and latitude axes are given as degrees . This should prevent a COARDS compliant application from mistaking the variables rlon und rlat to be actual longitude and latitude coordinates. The entries for these names in the standard name table indicate the appropriate sign conventions for the units of degrees .

An application can determine that x and y are the projection coordinates by recognizing the standard names projection_x_coordinate and projection_y_coordinate . The grid mapping variable Lambert_Conformal contains the mapping parameters as attributes, and is associated with the Temperature variable via its grid_mapping attribute.

5.6.1. Use of the CRS Well-known Text Format

An optional grid mapping attribute called crs_wkt may be used to specify multiple coordinate system properties in so-called well-known text format (usually abbreviated to CRS WKT or OGC WKT). The CRS WKT format is widely recognised and used within the geoscience software community. As such it represents a versatile mechanism for encoding information about a variety of coordinate reference system parameters in a highly compact notational form. The translation of CF coordinate variables to/from OGC Well-Known Text (WKT) format is shown in Examples 5.11 and 5.12 below and described in detail in https://github.com/cf-convention/cf-conventions/wiki/Mapping-from-CF-Grid-Mapping-Attributes-to-CRS-WKT-Elements.

Das crs_wkt attribute should comprise a text string that conforms to the WKT syntax as specified in reference [OGC_WKT-CRS]. If desired the text string may contain embedded newline characters to aid human readability. However, any such characters are purely cosmetic and do not alter the meaning of the attribute value. It is envisaged that the value of the crs_wkt attribute typically will be a single line of text, one intended primarily for machine processing. Other than the requirement to be a valid WKT string, the CF convention does not prescribe the content of the crs_wkt attribute since it will necessarily be context-dependent.

Wo ein crs_wkt attribute is added to a grid_mapping, the extended syntax for the grid_mapping attribute enables the list of variables containing coordinate values being referenced to be explicitly stated and the CRS WKT Axis order to be explicitly defined. The explicit definition of WKT CRS Axis order is expected by the OGC standards for referencing by coordinates. Software implementing these standards are likely to expect to receive coordinate value tuples, with the correct coordinate value order, along with the coordinate reference system definition that those coordinate values are defined with respect to.

The order of the <coordinatesVariable> references within the grid_mapping attribute definition defines the order of elements within a derived coordinate value tuple. This enables an application reading the data from a file to construct an array of coordinate value tuples, where each tuple is ordered to match the specification of the coordinate reference system being used whilst the array of tuples is structured according to the netCDF definition. It is the responsibility of the data producer to ensure that the <coordinatesVariable> list is consistent with the CRS WKT definition of CS AXIS, with the correct number of entries in the correct order (note: this is not a conformance requirement as CF conformance is not dependent on CRS WKT parsing).

For example, a file has two coordinate variables, lon and lat, and a grid mapping variable crs with an associated crs_wkt attribute the WKT definition defines the AXIS order as ["latitude", "longitude"]. The grid_mapping attribute is thus given a value crs:lat lon to define that where coordinate pairs are required, these shall be ordered (lat, lon), to be consistent with the provided crs_wkt string (and not order inverted). A 2-D array of (lat, lon) tuples can then be explicitly derived from the combination of the lat and lon variables.

Das crs_wkt attribute is intended to act as a Ergänzung to other single-property CF grid mapping attributes (as described in Appendix F) it is not intended to replace those attributes. If data producers omit the single-property grid mapping attributes in favour of the compound crs_wkt attribute, software which cannot interpret crs_wkt will be unable to use the grid_mapping information. Therefore the CRS should be described as thoroughly as possible with the single-property attributes as well as by crs_wkt .

There will be occasions when a given CRS property value is duplicated in both a single-property grid mapping attribute and the crs_wkt Attribut. In such cases the onus is on data producers to ensure that the property values are consistent. However, in those situations where two values of a given property are different, then the value specified by the single-property attribute shall take precedence. For example, if the semi-major axis length of the ellipsoid is defined by the grid mapping attribute semi_major_axis and also by the crs_wkt attribute (via the WKT SPHEROID[…​] element) then the former, being the more specific attribute, takes precedence. Naturally if the two values are equal then no ambiguity arises.

Likewise, in those cases where the value of a CRS WKT element should be used consistently across the CF-netCDF community (names of projections and projection parameters, for example) then, the values shown in https://github.com/cf-convention/cf-conventions/wiki/Mapping-from-CF-Grid-Mapping-Attributes-to-CRS-WKT-Elements should be preferred these are derived from the OGP/EPSG registry of geodetic parameters, which is considered to represent the definitive authority as regards CRS property names and values.

Examples 5.11 illustrates how the coordinate system properties specified via the crs grid mapping variable in Example 5.9 might be expressed using a crs_wkt Attribut. Example 5.12 also illustrates the addition of the crs_wkt attribute, but here the attribute is added to the crs variable of a simplified variant of Example 5.10. For brevity in Example 5.11, only the grid mapping variable and its grid_mapping_name and crs_wkt attributes are included all other elements are as per the Example 5.9. Names of projection PARAMETERs follow the spellings used in the EPSG geodetic parameter registry.

Example 5.12 illustrates how certain WKT elements - all of which are optional - can be used to specify CRS properties not covered by existing CF grid mapping attributes, including:

use of the VERT_DATUM element to specify vertical datum information

use of additional PARAMETER elements (albeit not essential ones in this example) to define the location of the false origin of the projection

use of AUTHORITY elements to specify object identifier codes assigned by an external authority, OGP/EPSG in this instance

Note: To enhance readability of these examples, the WKT value has been split across multiple lines and embedded quotation marks (") left unescaped - in real netCDF files such characters would need to be escaped. In CDL, within the CRS WKT definition string, newlines would need to be encoded within the string as and double quotes as " . Also for readability, we have dropped the quotation marks which would delimit the entire crs_wkt string. This pseudo CDL will not parse directly.

Note: There are unescaped double quotes and newlines and the quotation marks which would delimit the entire crs_wkt string are missing in this example. This is to enhance readability, but it means that this pseudo CDL will not parse directly.

The preceding two example (5.11 and 5.12) may be combined, if the data provider desires to provide explicit latitude and longitude coordinates as well as projection coordinates and to provide CRS WKT referencing for both sets of coordinates. This is demonstrated in example 5.13

Note: There are unescaped double quotes and newlines and the quotation marks which would delimit the entire crs_wkt string are missing in this example. This is to enhance readability, but it means that this pseudo CDL will not parse directly.

5.7. Scalar Coordinate Variables

When a variable has an associated coordinate which is single-valued, that coordinate may be represented as a scalar variable (i.e. a data variable which has no netCDF dimensions). Since there is no associated dimension these scalar coordinate variables should be attached to a data variable via the coordinates Attribut.

The use of scalar coordinate variables is a convenience feature which avoids adding size one dimensions to variables. A numeric scalar coordinate variable has the same information content and can be used in the same contexts as a size one numeric coordinate variable. Similarly, a string-valued scalar coordinate variable has the same meaning and purposes as a size one string-valued auxiliary coordinate variable (Section 6.1, "Labels"). Note however that use of this feature with a latitude, longitude, vertical, or time coordinate will inhibit COARDS conforming applications from recognizing them.

Once a name is used for a scalar coordinate variable it can not be used for a 1D coordinate variable. For this reason we strongly recommend against using a name for a scalar coordinate variable that matches the name of any dimension in the file.

If a data variable has two or more scalar coordinate variables, they are regarded as though they were all independent coordinate variables with dimensions of size one. If two or more single-valued coordinates are not independent, but have related values (this might be the case, for instance, for time and forecast period, or vertical coordinate and model level number, Section 6.2, "Alternative Coordinates"), they should be stored as coordinate or auxiliary coordinate variables of the same size one dimension, not as scalar coordinate variables.

In this example both the analysis time and the single pressure level are represented using scalar coordinate variables. The analysis time is identified by the standard name "forecast_reference_time" while the valid time of the forecast is identified by the standard name "time".


AreaDefinition Class Methods¶

There are four class methods available on the AreaDefinition class utilizing create_area_def() providing a simpler interface to the functionality described in the previous section. Hence each argument used below is the same as the create_area_def arguments described above and can be used in the same way (i.e. units). The following functions require area_id und projection along with a few other arguments:

From_extent¶

From_circle¶

From_area_of_interest¶

From_ul_corner¶


Beispiele

Display Regular Data Grid as Texture Map

Load elevation data and a geographic cells reference object for the Korean peninsula. Create a set of map axes for the Korean peninsula using worldmap .

Display the elevation data as a texture map. Apply a colormap appropriate for elevation data using demcmap .

Display Polygons on Map Projection

Import a shapefile containing the coastline coordinates of Africa, Asia, and Europe. Verify the data represents a polygon by querying its Geometry field.

Display the polygon on a world map. NaN values in coast separate the external continent boundary from the internal pond and lake boundaries.

Define Face Colors and Set Default Face Colors

Load sample data representing the USA. Set up an empty map axes with projection and limits suitable for displaying all 50 states.

Create a symbolization specification that sets the color of Alaska and Hawaii polygons to red.

Display all the other states, setting the default face color to blue and the default edge color to black.

Create Map and Display NaNs as Transparent

Load elevation data and a geographic cells reference object for the Korean peninsula. Insert a band of null values into the elevation data.

Create a set of map axes for the Korean peninsula using worldmap . Then, display the elevation data as a surface with transparent null values.

Display EGM96 Geoid Heights Masking Out Land Areas

Get geoid heights and a geographic postings reference object from the EGM96 geoid model. Then, display the geoid heights as a surface using an Eckert projection. Ensure the surface appears under the land mask by setting the 'CData' name-value pair to the geoid height data and the 'ZData' name-value pair to a matrix of zeros. Display the frame and grid of the map using framem and gridm .

Create a colorbar and add a text description. Then, mask out all the land.

Display EGM96 Geoid Heights as 3-D Surface

Get geoid heights and a geographic postings reference object from the EGM96 geoid model. Then, display the geoid heights as a surface using an Eckert projection.

Add light and material. Then, view the map as a 3-D surface.

Display Moon Albedo Using Orthographic Projection

Load moon albedo data and a geographic cells reference object.

Then, display the data. To do this, create a map axes object and specify its projection as orthographic. Display the data in the map axes as a texture map using the geoshow function. Then, change the colormap to grayscale and remove the axis lines.


NCO netCDF Operators

I have a NetCDF file structure as below. How can I delete the dimensions & variables such as "Band", "OpticalDepth", "ParticleType"?

I tried "ncks -x -v Band in.nc out.nc". It went thru from the command without any error. but when I ncdump out.nc, I see the variable/dimension "Band" is still in the file. Was ist los? Any suggestion how to do it?

dimensions:
Band = 1
time = UNLIMITED // (5 currently)
OpticalDepth = 1
ParticleType = 1
variables:
int Band(Band)
Band:long_name = "Band(fake)"
Band:units = "level"
Band:standard_name = "band"
double my_plottable_variable(time, OpticalDepth, Band, ParticleType)
my_plottable_variable:_FillValue = -9999.
my_plottable_variable:coordinates = "time"
my_plottable_variable:long_name = "my plottable variable"
int OpticalDepth(OpticalDepth)
OpticalDepth:long_name = "OpticalDepth(fake)"
OpticalDepth:units = "level"
OpticalDepth:standard_name = "opticaldepth"
int ParticleType(ParticleType)
ParticleType:long_name = "ParticleType(fake)"
ParticleType:units = "level"
ParticleType:standard_name = "particletype"
int time(time)
time:standard_name = "time"
time:units = "seconds since 1970-01-01 00:00:00"

Hoppla. just tried again, it complains with my command:
$ ncks -O -x -v new_misr.nc mytest.nc
ncks: ERROR mytest.nc neither exists locally nor matches remote filename patterns
$

so , still the same question -- what is the command to delete a dimension and a variable at the same time or separately?

ncks -C -O -x -v Band new_misr.nc mytest.nc

this will not delete the dimension because another variable (my_plottable_variable) needs it.

ncks -O -x -v Band,my_plottable_variable new_misr.nc mytest.nc

do appreciate. I see the problem now.

however, the variables 'my_plottable_variable' and "time" are to be used for a time series plot. so, deleting 'my_plottable_variable" is not a solution in my case. how to remove the "Band" dimension referenced by "my_plottable_variable"?

my final goal is to generate a file like below:

dimensions:
time = UNLIMITED // (5 currently)

variables:
double my_plottable_variable(time)
my_plottable_variable:_FillValue = -9999.
my_plottable_variable:coordinates = "time"
my_plottable_variable:long_name = "my plottable variable"

int time(time)
time:standard_name = "time"
time:units = "seconds since 1970-01-01 00:00:00"