[Python] Zeichnen Sie mit Plotly Höhendaten auf eine sphärische Oberfläche und zeichnen Sie einen Globus, der rund und rund gedreht werden kann

Einführung

Es unterstützt verschiedene Sprachen, einschließlich Python, und Sie können Plotly verwenden, mit dem Sie relativ einfach interaktive Diagramme zeichnen und Höhendaten verwenden können, um sich wie bei Google Earth zu drehen. Lassen Sie uns einen Globus erstellen.

Das Ergebnis sieht so aus. https://rkiuchir.github.io/3DSphericalTopo

3 Punkte

1. Erfassung und Lesen von Höhendaten
2. Konvertieren Sie die im orthogonalen Koordinatensystem ausgedrückten Breiten- und Längengradinformationen in das sphärische Koordinatensystem (Plotly Chart Studio: Heatmap-Plot auf einer sphärischen Karte. -siehe Handlung auf einer sphärischen Karte / # /))
3. Zeichnen Sie mit Plotly

Ausführungsumgebung

• macOS 10.14 Mojave
• Python 3.6 (Anaconda)
• Jupyter notebook

0. Plotly installieren

Es kann mit "pip install plotly" installiert werden.

1-1. Erfassung von Höhendaten

• Höhendaten (ETOPO1 Global Relief Model) mit einer Auflösung von 1 Minute weltweit sind in NOAA (National Oceanic and Atmospheric Administration) verfügbar. Diese Seite](https://www.ngdc.noaa.gov/mgg/global/) ist verfügbar.
• Laden Sie diesmal die im Raster registrierte Datei im netCDF-Format (ETOPO1_Ice_g_gdal.grd) von ETOPO1 Ice Surface herunter.

1-2. Höhendaten lesen

Lesen Sie zunächst die drei Daten für Breite / Länge und Höhe aus der Datendatei im netCDF-Format.

• Die hier definierte Eingabe ist der Breitengrad (lat_area), der Längengrad (lon_area) und die Auflösung (Auflösung) des angegebenen Bereichs.
• Die Ausgabe ist der Breitengrad (lon), der Längengrad (lat) und die Auflösung (topo) des Maschentyps im angegebenen Bereich.
• Der Auflösungsteil ist im Format ° (Grad) angegeben, und die höchste Auflösungsspezifikation ist 0,0167 °.

In Bezug auf die Auflösung werden die Daten gemäß dem angegebenen Auflösungswert übersprungen.

import numpy as np
from netCDF4 import Dataset

def Etopo(lon_area, lat_area, resolution):
    ### Input
    # resolution: resolution of topography for both of longitude and latitude [deg]
    # (Original resolution is 0.0167 deg)
    # lon_area and lat_area: the region of the map which you want like [100, 130], [20, 25]
    ###
    

    ### Output
    # Mesh type longitude, latitude, and topography data
    ### 

    # Read NetCDF data
    data = Dataset("ETOPO1_Ice_g_gdal.grd", "r")
    
    # Get data
    lon_range = data.variables['x_range'][:]
    lat_range = data.variables['y_range'][:]
    topo_range = data.variables['z_range'][:]
    spacing = data.variables['spacing'][:]
    dimension = data.variables['dimension'][:]
    z = data.variables['z'][:]

    lon_num = dimension[0]
    lat_num = dimension[1]
    
    # Prepare array
    lon_input = np.zeros(lon_num); lat_input = np.zeros(lat_num)
    for i in range(lon_num):
        lon_input[i] = lon_range[0] + i * spacing[0]

    for i in range(lat_num):
        lat_input[i] = lat_range[0] + i * spacing[1]

    # Create 2D array
    lon, lat = np.meshgrid(lon_input, lat_input)

    # Convert 2D array from 1D array for z value
    topo = np.reshape(z, (lat_num, lon_num))

    # Skip the data for resolution
    if ((resolution < spacing[0]) | (resolution < spacing[1])):
        print('Set the highest resolution')
    else:
        skip = int(resolution/spacing[0])
        lon = lon[::skip,::skip]
        lat = lat[::skip,::skip]
        topo = topo[::skip,::skip]

    topo = topo[::-1]

    # Select the range of map
    range1 = np.where((lon>=lon_area[0]) & (lon<=lon_area[1]))
    lon = lon[range1]; lat = lat[range1]; topo = topo[range1]
    range2 = np.where((lat>=lat_area[0]) & (lat<=lat_area[1]))
    lon = lon[range2]; lat = lat[range2]; topo = topo[range2]
    
    # Convert 2D again
    lon_num = len(np.unique(lon))
    lat_num = len(np.unique(lat))
    lon = np.reshape(lon, (lat_num, lon_num))
    lat = np.reshape(lat, (lat_num, lon_num))
    topo = np.reshape(topo, (lat_num, lon_num))
   
    return lon, lat, topo

2. Konvertieren Sie Längen- und Breitengradinformationen in sphärische Koordinaten

(Siehe Plotly Chart Studio: Heatmap-Plot auf einer sphärischen Karte)

Hier werden die Breiten- und Längengradinformationen, die durch die oben hergestellten orthogonalen Systemkoordinaten ausgedrückt werden, in das sphärische Koordinatensystem umgewandelt.

def degree2radians(degree):
    # convert degrees to radians
    return degree*np.pi/180

def mapping_map_to_sphere(lon, lat, radius=1):
    #this function maps the points of coords (lon, lat) to points onto the  
    sphere of radius radius
    
    lon=np.array(lon, dtype=np.float64)
    lat=np.array(lat, dtype=np.float64)
    lon=degree2radians(lon)
    lat=degree2radians(lat)
    xs=radius*np.cos(lon)*np.cos(lat)
    ys=radius*np.sin(lon)*np.cos(lat)
    zs=radius*np.sin(lat)
    return xs, ys, zs

3. Zeichnen Sie mit Plotly

Zeichnen wir nun die dreidimensionalen Daten von Breite, Länge und Höhe, die durch das sphärische Koordinatensystem in Plotly dargestellt werden.

Rufen Sie zunächst die in 1-2 vorbereitete Funktion auf, um die globalen Höhendaten zu lesen. Wenn Sie es mit einer zu hohen Auflösung lesen, erhöht sich die Datenmenge in der Größenordnung des Würfels. Diesmal wird die Auflösung auf 0,8 ° eingestellt.

# Import topography data
# Select the area you want
resolution = 0.8
lon_area = [-180., 180.]
lat_area = [-90., 90.]

# Get mesh-shape topography data
lon_topo, lat_topo, topo = ReadGeo.Etopo(lon_area, lat_area, resolution)

Konvertieren Sie es anschließend in ein sphärisches Koordinatensystem mit der in 2 vorbereiteten Funktion.

xs, ys, zs = mapping_map_to_sphere(lon_topo, lat_topo)

Und von hier aus werden wir tatsächlich mit dem Zeichnen fortfahren.

Definieren Sie zunächst die Farbskala, mit der die Höhendaten gezeichnet werden.

# Import color scale
import Plotly_code as Pcode
name = "topo"
Ctopo = Pcode.Colorscale_Plotly(name)
cmin = -8000
cmax = 8000

Zeichnen Sie dann mit Plotly. Hier geben Sie die Eingabedaten und die Farbskala ein.

topo_sphere=dict(type='surface',
            x=xs, 
            y=ys, 
            z=zs,
            colorscale=Ctopo,
            surfacecolor=topo,
            cmin=cmin, 
            cmax=cmax)
            )

Löschen Sie die Welle usw., damit sie gut aussieht.

noaxis=dict(showbackground=False,
            showgrid=False,
            showline=False,
            showticklabels=False,
            ticks='',
            title='',
            zeroline=False)

Verwenden Sie schließlich das Layout, um den Titel und die Hintergrundfarbe anzugeben. Diesmal ist die Hintergrundfarbe schwarz, mit ein wenig Bewusstsein für Google Earth.

import plotly.graph_objs as go

titlecolor = 'white'
bgcolor = 'black'

layout = go.Layout(
      autosize=False, width=1200, height=800,
      title = '3D spherical topography map',
      titlefont = dict(family='Courier New', color=titlecolor),
      showlegend = False,
      scene = dict(
        xaxis = noaxis,
        yaxis = noaxis,
        zaxis = noaxis,
        aspectmode='manual',
        aspectratio=go.layout.scene.Aspectratio(
            x=1, y=1, z=1)),
        paper_bgcolor = bgcolor,
        plot_bgcolor = bgcolor)

Zeichnen Sie dann mit der vorbereiteten (hier HTML-Ausgabe).

from plotly.offline import plot

plot_data=[topo_sphere]

fig = go.Figure(data=plot_data, layout=layout)
plot(fig, validate = False, filename='3DSphericalTopography.html', 
     auto_open=True)

abschließend

Ich denke, ich könnte damit eine Handlung eines Globus zeichnen, der sich wie am Anfang drehen lässt. Ich benutze es, indem ich die Verteilung der Erdbeben darüber lege.

Dieses Beispiel kann leicht angewendet werden, wenn auf einer Kugel gezeichnet wird, und ich denke, es hat eine breite Palette von Verwendungsmöglichkeiten. Darüber hinaus ist es nicht nur für Höhendaten nützlich, sondern auch zum Erstellen einer zweidimensionalen Karte, und es ist auch möglich, dreidimensionale Bilder zu zeichnen, die die Höhe mit Höhendaten ausdrücken.