main simulator file, incl. all sub-models and netCDF-interface
This commit is contained in:
parent
7fbb2b72f6
commit
aec1496620
379
main.py
379
main.py
|
|
@ -1,56 +1,274 @@
|
|||
import numpy as np
|
||||
import math
|
||||
import xarray as xr
|
||||
import matplotlib.pyplot as plt
|
||||
from scipy.spatial import cKDTree
|
||||
from scipy import interpolate
|
||||
import cartopy.crs as ccrs
|
||||
from datetime import datetime
|
||||
start = datetime.now()
|
||||
print(start)
|
||||
begin_time = datetime.now()
|
||||
from models.sun import sun_angles_analytical
|
||||
from models.sun import sun_angles_astropy
|
||||
from models.drag import drag, c_d
|
||||
from models.simple_atmosphere import T_air, p_air, rho_air
|
||||
from models.drag import drag #, c_d
|
||||
from input.user_input import *
|
||||
from input.natural_constants import *
|
||||
from astropy.time import Time
|
||||
import astropy.units as u
|
||||
import astropy.units as unit
|
||||
from models.thermal import AirMass
|
||||
from netCDF4 import Dataset
|
||||
import pandas as pd
|
||||
|
||||
data = pd.read_excel(r'C:\Users\marcel\PycharmProjects\MasterThesis\Mappe1.xls', sheet_name='Tabelle3')
|
||||
|
||||
comp_time = pd.DataFrame(data, columns= ['Time'])
|
||||
comp_height = pd.DataFrame(data, columns= ['Height'])
|
||||
comp_lat = pd.DataFrame(data, columns= ['Latitude'])
|
||||
comp_lon = pd.DataFrame(data, columns= ['Longitude'])
|
||||
|
||||
Reynolds = [0.0, 0.002, 0.004, 0.007, 0.01] #, 0.04, 0.07, 0.1, 0.2, 0.3, 0.5, 0.7, 1, 2, 3, 5, 7, 10, 20, 30, 50, 70, 100, 200,
|
||||
# 300, 500, 700, 1000, 2000, 3000, 5000, 7000, 10000, 20000, 30000, 50000]
|
||||
drag_coeff = [20000, 12000, 6000, 3429, 2400] #, 600, 343, 240, 120, 80, 49, 36.5, 26.5, 14.4, 10.4, 6.9, 5.4, 4.1, 2.55, 2, 1.5,
|
||||
# 1.27, 1.07, 0.77, 0.65, 0.55, 0.5, 0.46, 0.42, 0.4, 0.385, 0.39, 0.405, 0.45, 0.47, 0.49]
|
||||
|
||||
cd_interp1d = interpolate.interp1d(Reynolds, drag_coeff)
|
||||
|
||||
Re_list = []
|
||||
|
||||
def cd_func(Re):
|
||||
w = np.log10(Re)
|
||||
|
||||
if Re <= 0.01:
|
||||
c_d = 9/2 + 24/Re
|
||||
elif Re <= 20:
|
||||
c_d = 24/Re * (1 + 0.1315 * Re ** (0.82 - 0.05 * w))
|
||||
elif Re <= 260:
|
||||
c_d = 24/Re * (1 + 0.1935 * Re ** 0.6305)
|
||||
elif Re <= 1.5 * 10e3:
|
||||
c_d = 10 ** (1.6435 - 1.1242 * w + 0.1558 * w ** 2)
|
||||
elif Re <= 1.2 * 10e4:
|
||||
c_d = 10 ** (-2.4571 + 2.5558 * w - 0.9295 * w ** 2 + 0.1049 * w ** 3)
|
||||
elif Re <= 4.4 * 10e4:
|
||||
c_d = 10 ** (-1.9181 + 0.6370 * w - 0.0636 * w ** 3)
|
||||
elif Re <= 3.38 * 10e5:
|
||||
c_d = 10 ** (-4.3390 + 1.5809 * w - 0.1546 * w ** 2)
|
||||
elif Re <= 4 * 10e5:
|
||||
c_d = 29.78 - 5.3 * w
|
||||
elif Re <= 10 ** 6:
|
||||
c_d = 0.1 * w - 0.49
|
||||
elif Re > 10e6:
|
||||
c_d = 0.19 - 8 * 10e4/Re
|
||||
else:
|
||||
c_d = 0.47
|
||||
|
||||
return c_d
|
||||
|
||||
def cd_func_garg(Re):
|
||||
if Re == 0:
|
||||
c_d = 2.0
|
||||
elif Re < 2 * 10e5:
|
||||
c_d = 5.4856 * 10e9 * np.tanh(4.3774 * 10e-9 / Re) + 0.0709 * np.tanh(700.6575 / Re) + 0.3894 * np.tanh(74.1539 / Re) - 0.1198 * np.tanh(7429.0843 / Re) + 1.7174 * np.tanh(9.9851 / Re + 2.3384) + 0.4744
|
||||
elif 2 * 10e5 < Re <= 10e6:
|
||||
c_d = 8 * 10e-6 * ((Re/6530) ** 2 + np.tanh(Re) - 8 * np.log(Re) / np.log(10)) - 0.4119 * np.exp(-2.08 * 10 ** 43 / (Re + Re ** 2) ** 4) - 2.1344 * np.exp(-((np.log(Re ** 2 + 10.7563) / np.log(10)) ** 2 + 9.9867) / Re) + 0.1357 * np.exp(-((Re / 1620) ** 2 + 10370) / Re) - 8.5 * 10e-3 * (2 * np.log(np.tanh(np.tanh(Re))) / np.log(10) - 2825.7162) / Re + 2.4795
|
||||
elif Re > 10e6:
|
||||
c_d = 0.212546
|
||||
else:
|
||||
c_d = 2.0
|
||||
return 0.80 # c_d
|
||||
|
||||
#cd_interp1d = interpolate.interp1d(Reynolds, drag_coeff)
|
||||
|
||||
def transform(lon, lat, t):
|
||||
# WGS 84 reference coordinate system parameters
|
||||
A = 6378137.0 # major axis [m]
|
||||
E2 = 6.69437999014e-3 # eccentricity squared
|
||||
|
||||
t_s = 1000 * t
|
||||
|
||||
lon_rad = np.radians(lon)
|
||||
lat_rad = np.radians(lat)
|
||||
# convert to cartesian coordinates
|
||||
r_n = A / (np.sqrt(1 - E2 * (np.sin(lat_rad) ** 2)))
|
||||
x = r_n * np.cos(lat_rad) * np.cos(lon_rad)
|
||||
y = r_n * np.cos(lat_rad) * np.sin(lon_rad)
|
||||
z = r_n * (1 - E2) * np.sin(lat_rad)
|
||||
return x, y, z, t_s
|
||||
|
||||
data = Dataset("test2021.nc", "r", format="NETCDF4")
|
||||
|
||||
ERAtime = data.variables['time'][:] # time
|
||||
ERAlat = data.variables['latitude'][:] # latitude
|
||||
ERAlon = data.variables['longitude'][:] # longitude
|
||||
ERAz = data.variables['z'][:]/g # geopotential to geopotential height
|
||||
ERApress = data.variables['level'][:] # pressure level
|
||||
ERAtemp = data.variables['t'][:] # temperature in K
|
||||
|
||||
vw_x = data.variables['u'][:] # v_x
|
||||
vw_y = data.variables['v'][:] # v_y
|
||||
vw_z = data.variables['w'][:] # v_z
|
||||
|
||||
lon_era2d, lat_era2d, time_era = np.meshgrid(ERAlon, ERAlat, ERAtime)
|
||||
|
||||
xs, ys, zs, ts = transform(lon_era2d.flatten(), lat_era2d.flatten(), time_era.flatten())
|
||||
|
||||
tree = cKDTree(np.column_stack((xs, ys, zs, ts))) # !
|
||||
|
||||
|
||||
|
||||
# t: simulation time (0 .. t_max) in [s]
|
||||
# utc = Time(start_utc) + t * u.second # current time (UTC)
|
||||
|
||||
lat = start_lat # deg
|
||||
lon = start_lon # deg
|
||||
|
||||
date = Time('2020-12-13 12:00:00.000')
|
||||
|
||||
|
||||
# INITIALISATION
|
||||
|
||||
t = 0
|
||||
t = 0 # simulation time in seconds
|
||||
h = start_height
|
||||
utc = Time(start_utc)
|
||||
epoch_diff = (Time(start_utc).jd - Time('1900-01-01 00:00:00.0').jd) * 24.000000
|
||||
t_epoch = epoch_diff
|
||||
|
||||
|
||||
|
||||
xt, yt, zt, tt = transform(lon, lat, t) # test coordinates
|
||||
|
||||
d, inds = tree.query(np.column_stack((xt, yt, zt, tt)), k=8) # longitude, latitude, time in h # !
|
||||
w = 1.0 / d[0] ** 2
|
||||
|
||||
lat_sel = np.unravel_index(inds[0], lon_era2d.shape)[0]
|
||||
lon_sel = np.unravel_index(inds[0], lon_era2d.shape)[1]
|
||||
time_sel = np.unravel_index(inds[0], lon_era2d.shape)[2]
|
||||
|
||||
interp4d_height = np.ma.dot(w, ERAz[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
interp4d_temp = np.ma.dot(w, ERAtemp[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
interp4d_vw_x = np.ma.dot(w, vw_x[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
interp4d_vw_y = np.ma.dot(w, vw_y[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
interp4d_vw_z = np.ma.dot(w, vw_z[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
|
||||
pressure_hPa = np.array([1, 2, 3, 5, 7, 10, 20, 30, 50, 70, 100, 125, 150, 175, 200,225, 250, 300, 350, 400,
|
||||
450, 500, 550, 600, 650, 700, 750, 775, 800, 825, 850, 875, 900, 925, 950, 975, 1000])
|
||||
|
||||
pressure = 100 * pressure_hPa
|
||||
# print(pressure)
|
||||
|
||||
|
||||
# height_interp1d = interpolate.interp1d(interp4d_height, pressure)
|
||||
temp_interp1d = interpolate.interp1d(interp4d_height, interp4d_temp)
|
||||
press_interp1d = interpolate.interp1d(interp4d_height, pressure)
|
||||
vw_x_interp1d = interpolate.interp1d(interp4d_height, interp4d_vw_x)
|
||||
vw_y_interp1d = interpolate.interp1d(interp4d_height, interp4d_vw_y)
|
||||
vw_z_interp1d = interpolate.interp1d(interp4d_height, interp4d_vw_z)
|
||||
|
||||
|
||||
#u = 0
|
||||
#v = 0
|
||||
#w = 0
|
||||
T_air = temp_interp1d(h)
|
||||
p_air = press_interp1d(h)
|
||||
rho_air = p_air/(R_air * T_air)
|
||||
|
||||
u = vw_x_interp1d(h)
|
||||
v = vw_y_interp1d(h)
|
||||
w = -1 / g * vw_z_interp1d(h) * R_air * T_air / p_air
|
||||
|
||||
|
||||
v_x = 0
|
||||
v_y = 0
|
||||
v_z = 0
|
||||
p_gas = p_air(h)
|
||||
T_gas = T_air(h)
|
||||
|
||||
v_rel = 0
|
||||
|
||||
T_gas = T_air
|
||||
T_film = T_gas
|
||||
|
||||
t_list = []
|
||||
h_list = []
|
||||
v_list = []
|
||||
lat_list = []
|
||||
lon_list = []
|
||||
p_list = []
|
||||
Temp_list = []
|
||||
rho_list = []
|
||||
|
||||
while t <= t_end and h >= 0:
|
||||
t_list.append(t)
|
||||
h_list.append(h)
|
||||
v_list.append(v_rel)
|
||||
lat_list.append(lat)
|
||||
lon_list.append(lon)
|
||||
|
||||
p_list.append(p_air)
|
||||
Temp_list.append(T_air)
|
||||
rho_list.append(rho_air)
|
||||
|
||||
t_epoch = epoch_diff + t/3600
|
||||
|
||||
xt, yt, zt, tt = transform(lon, lat, t_epoch) # current balloon coordinates in cartesian coordinates: x [m], y [m], z [m], time [h since 1900-01-01 00:00:00.0]
|
||||
|
||||
d, inds = tree.query(np.column_stack((xt, yt, zt, tt)), k=8) # longitude, latitude, time in h # !
|
||||
w = 1.0 / d[0] ** 2
|
||||
|
||||
lat_sel = np.unravel_index(inds[0], lon_era2d.shape)[0]
|
||||
lon_sel = np.unravel_index(inds[0], lon_era2d.shape)[1]
|
||||
time_sel = np.unravel_index(inds[0], lon_era2d.shape)[2]
|
||||
|
||||
interp4d_height = np.ma.dot(w, ERAz[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
interp4d_temp = np.ma.dot(w, ERAtemp[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
interp4d_vw_x = np.ma.dot(w, vw_x[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
interp4d_vw_y = np.ma.dot(w, vw_y[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
interp4d_vw_z = np.ma.dot(w, vw_z[time_sel, :, lat_sel, lon_sel]) / np.sum(w)
|
||||
|
||||
press_interp1d = interpolate.interp1d(interp4d_height, pressure)
|
||||
temp_interp1d = interpolate.interp1d(interp4d_height, interp4d_temp)
|
||||
vw_x_interp1d = interpolate.interp1d(interp4d_height, interp4d_vw_x)
|
||||
vw_y_interp1d = interpolate.interp1d(interp4d_height, interp4d_vw_y)
|
||||
vw_z_interp1d = interpolate.interp1d(interp4d_height, interp4d_vw_z)
|
||||
|
||||
|
||||
try:
|
||||
p_air = press_interp1d(h)
|
||||
T_air = temp_interp1d(h)
|
||||
except:
|
||||
if (h > interp4d_height[0]):
|
||||
h = interp4d_height[0]
|
||||
else:
|
||||
h = interp4d_height[-1]
|
||||
|
||||
p_air = press_interp1d(h)
|
||||
T_air = temp_interp1d(h)
|
||||
|
||||
#print("height")
|
||||
#print(h)
|
||||
#print("temperature")
|
||||
#print(temp_interp1d(h))
|
||||
#print("pressure")
|
||||
#print(press_interp1d(h))
|
||||
#print("vw_x")
|
||||
#print(vw_x_interp1d(h))
|
||||
|
||||
p_gas = p_air
|
||||
rho_air = p_air / (R_air * T_air)
|
||||
u = vw_x_interp1d(h)
|
||||
v = vw_y_interp1d(h)
|
||||
w = -1/g * vw_z_interp1d(h) * R_air * T_air / p_air
|
||||
|
||||
v_relx = u - v_x
|
||||
v_rely = v - v_y
|
||||
v_relz = w - v_z
|
||||
|
||||
v_rel = (v_relx ** 2 + v_rely ** 2 + v_relz ** 2) ** (1/2)
|
||||
|
||||
rho_gas = p_gas/(R_gas * T_gas) # calculate gas density through ideal(!) gas equation
|
||||
|
||||
V_b = m_gas/rho_gas # calculate balloon volume from current gas mass and gas density
|
||||
|
||||
if V_b > V_design:
|
||||
V_b = V_design
|
||||
m_gas = V_design * rho_gas
|
||||
V_b = V_design
|
||||
m_gas = V_design * rho_gas
|
||||
else:
|
||||
pass
|
||||
|
||||
m_gross = m_pl + m_film
|
||||
m_tot = m_pl + m_film + m_gas
|
||||
m_virt = m_tot + c_virt * rho_air(h) * V_b
|
||||
m_virt = m_tot + c_virt * rho_air * V_b
|
||||
|
||||
d_b = 1.383 * V_b ** (1/3) # calculate diameter of balloon from its volume
|
||||
L_goreB = 1.914 * V_b ** (1/3)
|
||||
|
|
@ -60,13 +278,6 @@ while t <= t_end and h >= 0:
|
|||
A_eff = 0.65 * A_surf + 0.35 * A_surf1
|
||||
A_top = np.pi/4 * d_b ** 2
|
||||
|
||||
D = drag(c_d, rho_air(h), d_b, v_z) # calculate drag force
|
||||
I = g * V_b * (rho_air(h) - rho_gas) # calculate gross inflation
|
||||
W = g * m_gross # calculate weight (force)
|
||||
F = I - W - D * np.sign(v_z)
|
||||
|
||||
v_zrel = v_z
|
||||
|
||||
AZ, ELV = sun_angles_analytical(lat, lon, utc)
|
||||
|
||||
A_proj = A_top * (0.9125 + 0.0875 * np.cos(np.pi - 2 * np.deg2rad(ELV)))
|
||||
|
|
@ -74,8 +285,8 @@ while t <= t_end and h >= 0:
|
|||
# CALCULATIONS FOR THERMAL MODEL
|
||||
|
||||
if ELV >= -(180 / np.pi * np.arccos(R_E / (R_E + h))):
|
||||
tau_atm = 0.5 * (np.exp(-0.65 * AirMass(p_air(h), p_0, ELV, h)) + np.exp(-0.095 * AirMass(p_air(h), p_0, ELV, h)))
|
||||
tau_atmIR = 1.716 - 0.5 * (np.exp(-0.65 * p_air(h) / p_0) + np.exp(-0.095 * p_air(h) / p_0))
|
||||
tau_atm = 0.5 * (np.exp(-0.65 * AirMass(p_air, p_0, ELV, h)) + np.exp(-0.095 * AirMass(p_air, p_0, ELV, h)))
|
||||
tau_atmIR = 1.716 - 0.5 * (np.exp(-0.65 * p_air / p_0) + np.exp(-0.095 * p_air / p_0))
|
||||
else:
|
||||
tau_atm = 0
|
||||
tau_atmIR = 0
|
||||
|
|
@ -95,17 +306,19 @@ while t <= t_end and h >= 0:
|
|||
else:
|
||||
q_Albedo = Albedo * I_Sun * np.sin(np.deg2rad(ELV))
|
||||
|
||||
my_air = (1.458 * 10 ** -6 * T_air(h) ** 1.5) / (T_air(h) + 110.4)
|
||||
my_air = (1.458 * 10 ** -6 * T_air ** 1.5) / (T_air + 110.4)
|
||||
my_gas = 1.895 * 10 ** -5 * (T_gas / 273.15) ** 0.647
|
||||
k_air = 0.0241 * (T_air(h) / 273.15) ** 0.9
|
||||
k_air = 0.0241 * (T_air / 273.15) ** 0.9
|
||||
k_gas = 0.144 * (T_gas / 273.15) ** 0.7
|
||||
Pr_air = 0.804 - 3.25 * 10 ** (-4) * T_air(h)
|
||||
Pr_air = 0.804 - 3.25 * 10 ** (-4) * T_air
|
||||
Pr_gas = 0.729 - 1.6 * 10 ** (-4) * T_gas
|
||||
|
||||
Gr_air = (rho_air(h) ** 2 * g * np.abs(T_film - T_air(h)) * d_b ** 3) / (T_air(h) * my_air ** 2)
|
||||
Gr_air = (rho_air ** 2 * g * np.abs(T_film - T_air) * d_b ** 3) / (T_air * my_air ** 2)
|
||||
Nu_air = 2 + 0.45 * (Gr_air * Pr_air) ** 0.25
|
||||
HC_free = Nu_air * k_air / d_b
|
||||
Re = np.abs(v_zrel) * d_b * rho_air(h) / my_air
|
||||
Re = np.abs(v_relz) * d_b * rho_air / my_air
|
||||
|
||||
Re_list.append(Re)
|
||||
|
||||
HC_forced = k_air / d_b * (2 + 0.41 * Re ** 0.55)
|
||||
HC_internal = 0.13 * k_gas * ((rho_gas ** 2 * g * np.abs(T_film - T_gas) * Pr_gas) / (T_gas * my_air ** 2)) ** (
|
||||
|
|
@ -120,39 +333,121 @@ while t <= t_end and h >= 0:
|
|||
Q_IREarth = alpha_IR * A_surf * q_IREarth * ViewFactor * (1 + tau_IR / (1 - r_IR))
|
||||
Q_IRfilm = sigma * epsilon * alpha_IR * A_surf * T_film ** 4 * 1 / (1 - r_IR)
|
||||
Q_IRout = sigma * epsilon * A_surf * T_film ** 4 * (1 * tau_IR / (1 - r_IR))
|
||||
Q_ConvExt = HC_external * A_eff * (T_air(h) - T_film)
|
||||
Q_ConvExt = HC_external * A_eff * (T_air - T_film)
|
||||
Q_ConvInt = HC_internal * A_eff * (T_film - T_gas)
|
||||
|
||||
RoC = -v_z
|
||||
try:
|
||||
c_d = cd_func(Re)
|
||||
except:
|
||||
print(Re)
|
||||
|
||||
D = drag(c_d, rho_air, d_b, v_rel) # calculate drag force
|
||||
|
||||
if v_rel == 0:
|
||||
Drag_x = 0
|
||||
Drag_y = 0
|
||||
Drag_z = 0
|
||||
else:
|
||||
Drag_x = D * v_relx / v_rel
|
||||
Drag_y = D * v_rely / v_rel
|
||||
Drag_z = D * v_relz / v_rel
|
||||
|
||||
I = g * V_b * (rho_air - rho_gas) # calculate gross inflation
|
||||
W = g * m_gross # calculate weight (force)
|
||||
F = I - W + Drag_z
|
||||
|
||||
# RoC = -v_z
|
||||
# dT_gas = (Q_ConvInt / (gamma * m_gas * c_v) - (gamma - 1) / gamma * rho_air(h) * g / (rho_gas * R_gas) * RoC) * dt
|
||||
# dT_film = ((Q_Sun + Q_Albedo + Q_IREarth + Q_IRfilm + Q_ConvExt - Q_ConvInt - Q_IRout) / (c_f * m_film)) * dt
|
||||
# v_w = w
|
||||
# v = v_z
|
||||
|
||||
v_w = 0
|
||||
v = v_z
|
||||
s = h
|
||||
|
||||
a = F/m_virt
|
||||
a_x = Drag_x/m_virt
|
||||
a_y = Drag_y/m_virt
|
||||
a_z = F/m_virt
|
||||
|
||||
s_n = s + dt * v + 0.5 * a * dt ** 2
|
||||
# DIFFERENTIAL EQUATIONS
|
||||
|
||||
dx = dt * v_x + 0.5 * a_x * dt ** 2 # lon
|
||||
dy = dt * v_y + 0.5 * a_y * dt ** 2 # lat
|
||||
|
||||
lon = lon + np.rad2deg(np.arctan(dx/((6371229.0 + h) * np.cos(np.deg2rad(lat)))))
|
||||
lat = lat + np.rad2deg(np.arctan(dy/(6371229.0 + h)))
|
||||
|
||||
v_xn = v_x + dt * a_x
|
||||
v_yn = v_y + dt * a_y
|
||||
|
||||
s_n = s + dt * v_z + 0.5 * a_z * dt ** 2
|
||||
dh = s_n - s
|
||||
v_n = v - v_w + dt * a
|
||||
v_n = v_z + dt * a_z
|
||||
|
||||
T_gn = T_gas + dt * (Q_ConvInt / (gamma * c_v * m_gas) - g * R_gas * T_gas * dh / (c_v * gamma * T_air(s) * R_air))
|
||||
T_gn = T_gas + dt * (Q_ConvInt / (gamma * c_v * m_gas) - g * R_gas * T_gas * dh / (c_v * gamma * T_air * R_air))
|
||||
T_en = T_film + (Q_Sun + Q_Albedo + Q_IREarth + Q_IRfilm + Q_ConvExt - Q_ConvInt - Q_IRout) / (c_f * m_film) * dt
|
||||
|
||||
T_g = T_gn
|
||||
T_e = T_en
|
||||
s = s_n
|
||||
v = v_n
|
||||
v_x = v_xn
|
||||
v_y = v_yn
|
||||
v_z = v_n
|
||||
|
||||
h = s
|
||||
T_gas = T_g
|
||||
T_film = T_e
|
||||
v_z = v
|
||||
# v_z = v
|
||||
|
||||
t += dt # time increment
|
||||
utc = Time(start_utc) + t * u.second
|
||||
utc = Time(start_utc) + t * unit.second
|
||||
|
||||
plt.plot(t_list, h_list)
|
||||
|
||||
|
||||
#print(len(lon_list), len(lat_list))
|
||||
#"""
|
||||
ax = plt.axes(projection=ccrs.Mollweide())
|
||||
#ax.coastlines()
|
||||
ax.stock_img()
|
||||
plt.plot(start_lon, start_lat, 'rx', transform=ccrs.Geodetic())
|
||||
plt.plot(lon_list, lat_list, 'k--', transform=ccrs.Geodetic())
|
||||
plt.plot(comp_lon, comp_lat, 'r-', transform=ccrs.Geodetic())
|
||||
#plt.xticks()
|
||||
plt.show()
|
||||
#"""
|
||||
|
||||
|
||||
### WORKS!
|
||||
###
|
||||
###dset = xr.open_dataset("test2021.nc")
|
||||
###
|
||||
###print(dset['t'][1][0]) # [time] [level]
|
||||
###
|
||||
#fig = plt.figure() #figsize=[120,50])
|
||||
###
|
||||
#ax = fig.add_subplot(111, projection=ccrs.PlateCarree())
|
||||
###
|
||||
###dset['t'][1][0].plot(ax=ax, cmap='jet', transform=ccrs.PlateCarree())
|
||||
###ax.coastlines()
|
||||
###plt.show()
|
||||
|
||||
#print(len(ERAlat))
|
||||
#print(len(ERAlon))
|
||||
#print(len(ERAz))
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
#"""
|
||||
plt.plot(comp_time, comp_height, 'r--', label='PoGo Flight Test')
|
||||
#plt.plot(t_list, Re_list, 'k-', label="Reynolds-Number")
|
||||
#plt.plot(t_list, v_list, 'k--', label='Ballon v_rel')
|
||||
plt.plot(t_list, h_list, 'k-', label='Balloon Altitude')
|
||||
#plt.plot(t_list, p_list, 'r-', label='Air Pressure [Pa]')
|
||||
#plt.plot(t_list, Temp_list, 'b-', label='Air Temperature [K]')
|
||||
#plt.plot(t_list, rho_list, 'g-', label='Air Density in [kg/m$^3$]')
|
||||
plt.xlabel('time in s')
|
||||
plt.ylabel('Balloon Altitude in m')
|
||||
plt.legend()
|
||||
plt.show()
|
||||
#"""
|
||||
|
|
|
|||
Loading…
Reference in New Issue