From 2f8bdb612022281101977062511929a4209b727b Mon Sep 17 00:00:00 2001 From: Harrison Nicholls Date: Thu, 18 Jul 2024 11:55:48 +0100 Subject: [PATCH 1/5] Merge sulfur fixes and plotting updates. Downloads stellar spectra from OSF --- README.md | 2 +- examples/SocRadConv.py | 4 +- examples/demo_instellation.py | 4 +- examples/demo_runaway_greenhouse.py | 4 +- src/janus/modules/plot_flux_balance.py | 168 ++----------------------- src/janus/utils/GeneralAdiabat.py | 26 ++-- src/janus/utils/data.py | 39 ++++-- src/janus/utils/phys.py | 4 +- tests/test_instellation.py | 4 +- tests/test_runaway_greenhouse.py | 4 +- 10 files changed, 66 insertions(+), 193 deletions(-) diff --git a/README.md b/README.md index d03a7a8..ac7065b 100644 --- a/README.md +++ b/README.md @@ -1,4 +1,4 @@ -## JANUS (temperature structure generator) +## JANUS (1D convective atmosphere model) Generates a temperature profile using the generalised moist pseudoadiabat and a prescribed stratosphere. Calculates radiative fluxes using SOCRATES. diff --git a/examples/SocRadConv.py b/examples/SocRadConv.py index e95446a..5b4cd40 100755 --- a/examples/SocRadConv.py +++ b/examples/SocRadConv.py @@ -21,7 +21,7 @@ from importlib.resources import files from janus.modules import RadConvEqm, plot_fluxes, plot_emission -from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, plot_adiabats, ReadBandEdges, StellarSpectrum +from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, DownloadStellarSpectra, plot_adiabats, ReadBandEdges, StellarSpectrum import mors #################################### @@ -70,7 +70,7 @@ #Download required spectral files DownloadSpectralFiles("/Dayspring") - DownloadSpectralFiles("/stellar_spectra") + DownloadStellarSpectra() # Move/prepare spectral file print("Inserting stellar spectrum") diff --git a/examples/demo_instellation.py b/examples/demo_instellation.py index 5b7c97b..160b3d1 100755 --- a/examples/demo_instellation.py +++ b/examples/demo_instellation.py @@ -12,7 +12,7 @@ import numpy as np from janus.modules import MCPA_CBL -from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, ReadBandEdges, StellarSpectrum +from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, DownloadStellarSpectra, ReadBandEdges, StellarSpectrum import mors @@ -37,7 +37,7 @@ #Download required spectral files DownloadSpectralFiles("/Oak") - DownloadSpectralFiles("/stellar_spectra") + DownloadStellarSpectra() # Setup spectral file print("Inserting stellar spectrum") diff --git a/examples/demo_runaway_greenhouse.py b/examples/demo_runaway_greenhouse.py index b870f5b..5b03165 100755 --- a/examples/demo_runaway_greenhouse.py +++ b/examples/demo_runaway_greenhouse.py @@ -11,7 +11,7 @@ from importlib.resources import files from janus.modules import RadConvEqm -from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, ReadBandEdges, StellarSpectrum +from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, DownloadStellarSpectra, ReadBandEdges, StellarSpectrum import mors if __name__=='__main__': @@ -36,7 +36,7 @@ #Download required spectral files DownloadSpectralFiles("/Oak") - DownloadSpectralFiles("/stellar_spectra") + DownloadStellarSpectra() # Setup spectral file print("Inserting stellar spectrum") diff --git a/src/janus/modules/plot_flux_balance.py b/src/janus/modules/plot_flux_balance.py index 4fb58e8..cdc7c60 100644 --- a/src/janus/modules/plot_flux_balance.py +++ b/src/janus/modules/plot_flux_balance.py @@ -16,155 +16,6 @@ from janus.modules.spectral_planck_surface import surf_Planck_nu import janus.utils.GeneralAdiabat as ga # Moist adiabat with multiple condensibles - -def plot_flux_balance(atm_dry, atm_moist, cp_dry, time, dirs): - - fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(13,10)) - # sns.set_style("ticks") - # sns.despine() - - # Line settings - col_idx = 3 - col_vol1 = "H2O" - col_vol2 = "N2" - col_vol3 = "H2" - col_vol4 = "O2" - - # Temperature vs. pressure - ax1.semilogy(atm_moist.tmp,atm_moist.p, color=ga.vol_colors[col_vol1][col_idx+1], ls="-", label=r'Moist adiabat') - if cp_dry == True: ax1.semilogy(atm_dry.tmp,atm_dry.p, color=ga.vol_colors[col_vol3][col_idx+1], ls="-", label=r'Dry adiabat') - ax1.legend() - ax1.invert_yaxis() - ax1.set_xlabel(r'Temperature $T$ (K)') - ax1.set_ylabel(r'Pressure $P$ (Pa)') - # ax1.set_ylim(bottom=atm_moist.ps*1.01) - ax1.set_ylim(top=atm_moist.ptop, bottom=atm_moist.ps) - - # Print active species - active_species = r"" - for vol in atm_moist.vol_list: - if atm_moist.vol_list[vol] > 1e-5: - active_species = active_species + ga.vol_latex[vol] + ", " - active_species = active_species[:-2] - ax1.text(0.02, 0.02, r"Active species: "+active_species, va="bottom", ha="left", fontsize=10, transform=ax1.transAxes, bbox=dict(fc='white', ec="white", alpha=0.5, boxstyle='round', pad=0.1), color=ga.vol_colors["black_1"] ) - - # Fluxes vs. pressure - - # Zero line - ax2.axvline(0, color=ga.vol_colors["qgray_light"], lw=0.5) - # ax2.axvline(-1e+3, color=ga.vol_colors["qgray_light"], lw=0.5) - # ax2.axvline(1e+3, color=ga.vol_colors["qgray_light"], lw=0.5) - - # SW down / instellation flux - if cp_dry == True: ax2.semilogy(atm_dry.SW_flux_down*(-1),atm_dry.pl, color=ga.vol_colors[col_vol3][col_idx], ls=":") - ax2.semilogy(atm_moist.SW_flux_down*(-1),atm_moist.pl, color=ga.vol_colors[col_vol2][col_idx], ls=":", label=r'$F_{\odot}^{\downarrow}$') - - # LW down / thermal downward flux - if cp_dry == True: ax2.semilogy(atm_dry.LW_flux_down*(-1),atm_dry.pl, color=ga.vol_colors[col_vol3][col_idx], ls="--") - ax2.semilogy(atm_moist.LW_flux_down*(-1),atm_moist.pl, color=ga.vol_colors[col_vol2][col_idx+1], ls="--", label=r'$F_\mathrm{t}^{\downarrow}$') - # ls=(0, (3, 1, 1, 1)) - - # Thermal upward flux, total - if cp_dry == True: ax2.semilogy(atm_dry.flux_up_total,atm_dry.pl, color=ga.vol_colors[col_vol3][col_idx], ls="--") - ax2.semilogy(atm_moist.flux_up_total,atm_moist.pl, color=ga.vol_colors[col_vol1][col_idx], ls="--", label=r'$F_\mathrm{t}^{\uparrow}$') - - # Net flux - if cp_dry == True: ax2.semilogy(atm_dry.net_flux,atm_dry.pl, color=ga.vol_colors[col_vol3][6], ls="-", lw=2) - ax2.semilogy(atm_moist.net_flux,atm_moist.pl, color=ga.vol_colors[col_vol1][6], ls="-", lw=2, label=r'$F_\mathrm{net}^{\uparrow}$') - - # # SW up - # if cp_dry == True: ax2.semilogy(atm_dry.SW_flux_up,atm_dry.pl, color=ga.vol_colors[col_vol3][col_idx], ls=":") - # ax2.semilogy(atm_moist.SW_flux_up,atm_moist.pl, color=ga.vol_colors[col_vol1][col_idx], ls=":", label=r'$F_\mathrm{SW}^{\uparrow}$') - - # # LW up - # if cp_dry == True: ax2.semilogy(atm_dry.LW_flux_up,atm_dry.pl, color=ga.vol_colors[col_vol3][col_idx], ls=(0, (5, 1))) - # ax2.semilogy(atm_moist.LW_flux_up,atm_moist.pl, color=ga.vol_colors[col_vol1][col_idx], ls=(0, (5, 1)), label=r'$F_\mathrm{LW}^{\uparrow}$') - - ax2.legend(ncol=6, fontsize=10, loc=3) - ax2.invert_yaxis() - ax2.set_xscale("symlog") # https://stackoverflow.com/questions/3305865/what-is-the-difference-between-log-and-symlog - ax2.set_xlabel(r'Outgoing flux $F^{\uparrow}$ (W m$^{-2}$)') - ax2.set_ylabel(r'Pressure $P$ (Pa)') - ax2.set_ylim(top=atm_moist.ptop, bottom=atm_moist.ps) - - # Wavenumber vs. OLR - ax3.plot(atm_moist.band_centres, surf_Planck_nu(atm_moist)/atm_moist.band_widths, color="gray",ls='--',label=str(round(atm_moist.ts))+' K blackbody') - if cp_dry == True: ax3.plot(atm_dry.band_centres, atm_dry.net_spectral_flux[:,0]/atm_dry.band_widths, color=ga.vol_colors[col_vol3][col_idx]) - ax3.plot(atm_moist.band_centres, atm_moist.net_spectral_flux[:,0]/atm_moist.band_widths, color=ga.vol_colors[col_vol1][col_idx+1]) - ax3.set_ylabel(r'Spectral flux density (W m$^{-2}$ cm$^{-1}$)') - ax3.set_xlabel(r'Wavenumber (cm$^{-1}$)') - ax3.legend(loc=1) - ymax_plot = 1.2*np.max(atm_moist.net_spectral_flux[:,0]/atm_moist.band_widths) - ax3.set_ylim(bottom=0, top=ymax_plot) - ax3.set_xlim(left=0, right=np.max(np.where(atm_moist.net_spectral_flux[:,0]/atm_moist.band_widths > ymax_plot/1000., atm_moist.band_centres, 0.))) - # ax3.set_xlim(left=0, right=30000) - - - # Heating versus pressure - ax4.axvline(0, color=ga.vol_colors["qgray_light"], lw=0.5) - - if cp_dry == True: - ax4.plot(atm_dry.LW_heating, atm_dry.p, ls="--", color=ga.vol_colors[col_vol3][col_idx+1]) - ax4.plot(atm_dry.net_heating, atm_dry.p, lw=2, color=ga.vol_colors[col_vol3][col_idx+1]) - ax4.plot(atm_dry.SW_heating, atm_dry.p, ls=":", color=ga.vol_colors[col_vol3][col_idx+1]) - - ax4.plot(atm_moist.LW_heating, atm_moist.p, ls="--", color=ga.vol_colors[col_vol1][col_idx+1], label=r'LW') - ax4.plot(atm_moist.net_heating, atm_moist.p, lw=2, color=ga.vol_colors[col_vol1][col_idx+1], label=r'Net') - ax4.plot(atm_moist.SW_heating, atm_moist.p, ls=":", color=ga.vol_colors[col_vol1][col_idx+1], label=r'SW') - - # Plot tropopause - trpp_idx = int(atm_moist.trppidx) - if trpp_idx > 0: - ax4.axhline(atm_moist.pl[trpp_idx], color=ga.vol_colors[col_vol2][col_idx], lw=1.0, ls="-.", label=r'Tropopause') - - ax4.invert_yaxis() - ax4.legend(ncol=4, fontsize=10, loc=3) - ax4.set_ylabel(r'Pressure $P$ (Pa)') - ax4.set_xlabel(r'Heating (K/day)') - # ax4.set_xscale("log") - ax4.set_yscale("log") - x_minmax = np.max([abs(np.min(atm_moist.net_heating[10:])), abs(np.max(atm_moist.net_heating[10:]))]) - x_minmax = np.max([ 20, x_minmax ]) - if not math.isnan(x_minmax): - ax4.set_xlim(left=-x_minmax, right=x_minmax) - ax4.set_ylim(top=atm_moist.ptop, bottom=atm_moist.ps) - # ax4.set_yticks([1e-10, 1e-5, 1e0, 1e5]) - # ax4.set_xticks([0.1, 0.3, 1, 3, 10, 30, 100]) - # ax4.set_xticklabels(["0.1", "0.3", "1", "3", "10", "30", "100"]) - - # # Wavelength versus OLR log plot - # OLR_cm_moist = atm_moist.LW_spectral_flux_up[:,0]/atm_moist.band_widths - # wavelength_moist = [ 1e+4/i for i in atm_moist.band_centres ] # microns - # OLR_micron_moist = [ 1e+4*i for i in OLR_cm_moist ] # microns - # if cp_dry == True: - # OLR_cm_dry = atm_dry.LW_spectral_flux_up[:,0]/atm_dry.band_widths - # wavelength_dry = [ 1e+4/i for i in atm_dry.band_centres ] # microns - # OLR_micron_dry = [ 1e+4*i for i in OLR_cm_dry ] # microns - # ax4.plot(wavelength_dry, OLR_micron_dry, color=ga.vol_colors[col_vol3][col_idx+1]) - - # ax4.plot(wavelength_moist, OLR_micron_moist, color=ga.vol_colors[col_vol1][col_idx+1]) - # ax4.set_ylabel(r'Spectral flux density (W m$^{-2}$ $\mu$m$^{-1}$)') - # ax4.set_xlabel(r'Wavelength $\lambda$ ($\mu$m)') - # ax4.set_xscale("log") - # ax4.set_yscale("log") - # ax4.set_xlim(left=0.1, right=100) - # ax4.set_ylim(bottom=1e-20, top=1e5) - # # ax4.set_yticks([1e-10, 1e-5, 1e0, 1e5]) - # ax4.set_xticks([0.1, 0.3, 1, 3, 10, 30, 100]) - # ax4.set_xticklabels(["0.1", "0.3", "1", "3", "10", "30", "100"]) - - fig.savefig(dirs["output"]+"/"+"TP_"+str(round(time["planet"]))+'.pdf', bbox_inches="tight") - plt.close(fig) - - # with open(dirs["output"]+"/"+str(int(time["planet"]))+"_atm.pkl", "wb") as atm_file: - # pkl.dump(atm, atm_file, protocol=pkl.HIGHEST_PROTOCOL) - - # # Save atm object to .json file - # json_atm = json.dumps(atm.__dict__) - # with open(dirs["output"]+"/"+str(int(time_current))+"_atm.json", "wb") as atm_file: - # json.dump(json_atm, atm_file) - - # Plotting def plot_fluxes(atm,filename='output/fluxes.pdf'): @@ -172,24 +23,23 @@ def plot_fluxes(atm,filename='output/fluxes.pdf'): ax.axvline(0,color='black',lw=0.8) - # for p in np.logspace(-6, 3, 80)*1.0e5: - # ax.axhline(y=p, color='purple', lw=1, alpha=0.3) + arr_p = atm.pl*1e-5 - ax.plot(atm.flux_up_total,atm.pl,color='red',label='UP',lw=1) - ax.plot(atm.SW_flux_up ,atm.pl,color='red',label='UP SW',linestyle='dotted',lw=2) - ax.plot(atm.LW_flux_up ,atm.pl,color='red',label='UP LW',linestyle='dashed',lw=1) + ax.plot(atm.flux_up_total,arr_p,color='red',label='UP',lw=1) + ax.plot(atm.SW_flux_up ,arr_p,color='red',label='UP SW',linestyle='dotted',lw=2) + ax.plot(atm.LW_flux_up ,arr_p,color='red',label='UP LW',linestyle='dashed',lw=1) - ax.plot(-1.0*atm.flux_down_total,atm.pl,color='green',label='DN',lw=2) - ax.plot(-1.0*atm.SW_flux_down ,atm.pl,color='green',label='DN SW',linestyle='dotted',lw=3) - ax.plot(-1.0*atm.LW_flux_down ,atm.pl,color='green',label='DN LW',linestyle='dashed',lw=2) + ax.plot(-1.0*atm.flux_down_total,arr_p,color='green',label='DN',lw=2) + ax.plot(-1.0*atm.SW_flux_down ,arr_p,color='green',label='DN SW',linestyle='dotted',lw=3) + ax.plot(-1.0*atm.LW_flux_down ,arr_p,color='green',label='DN LW',linestyle='dashed',lw=2) - ax.plot(atm.net_flux ,atm.pl,color='black',label='NET') + ax.plot(atm.net_flux ,arr_p,color='black',label='NET') ax.set_xlabel("Upward-directed flux [W m-2]") ax.set_xscale("symlog") ax.set_ylabel("Pressure [Pa]") ax.set_yscale("log") - ax.set_ylim([atm.pl[-1],atm.pl[0]]) + ax.set_ylim(top=np.amin(arr_p)/1.5, bottom=np.amax(arr_p)*1.5) ax.legend(loc='upper left') diff --git a/src/janus/utils/GeneralAdiabat.py b/src/janus/utils/GeneralAdiabat.py index 1b037e0..11fbec9 100644 --- a/src/janus/utils/GeneralAdiabat.py +++ b/src/janus/utils/GeneralAdiabat.py @@ -33,17 +33,19 @@ # https://chrisalbon.com/python/data_visualization/seaborn_color_palettes/ no_colors = 7 vol_colors = { - "H2O" : [mpl.colormaps["PuBu"](i) for i in np.linspace(0,1.0,no_colors)], - "CO2" : [mpl.colormaps["Reds"](i) for i in np.linspace(0,1.0,no_colors)], - "H2" : [mpl.colormaps["Greens"](i) for i in np.linspace(0,1.0,no_colors)], - "N2" : [mpl.colormaps["Purples"](i) for i in np.linspace(0,1.0,no_colors)], "O2" : [mpl.colormaps["Wistia"](i) for i in np.linspace(0,1.0,no_colors+2)], "O3" : [mpl.colormaps["spring"](i) for i in np.linspace(0,1.0,no_colors+2)], - "CH4" : [mpl.colormaps["RdPu"](i) for i in np.linspace(0,1.0,no_colors)], - "CO" : [mpl.colormaps["pink_r"](i) for i in np.linspace(0,1.0,no_colors)], "S" : [mpl.colormaps["YlOrBr"](i) for i in np.linspace(0,1.0,no_colors)], - "He" : [mpl.colormaps["Greys"](i) for i in np.linspace(0,1.0,no_colors)], - "NH3" : [mpl.colormaps["summer"](i) for i in np.linspace(0,1.0,no_colors)], + "H2O": "#027FB1", + "CO2": "#D24901", + "H2" : "#008C01", + "CH4": "#C720DD", + "CO" : "#D1AC02", + "N2" : "#870036", + "S2" : "#FF8FA1", + "SO2": "#00008B", + "He" : "#30FF71", + "NH3": "#675200", "mixtures" : [mpl.colormaps["Set3"](i) for i in np.linspace(0,1.0,no_colors)], "H2O-CO2" : mpl.colormaps["Set3"](1.0/no_colors), "CO2-H2O" : mpl.colormaps["Set3"](1.0/no_colors), @@ -943,17 +945,17 @@ def plot_adiabats(atm,filename='output/general_adiabat.pdf'): # Saturation vapor pressure for given temperature Psat_array = [ p_sat(vol, T,water_lookup=atm.water_lookup) for T in T_sat_array ] - ax1.semilogy( T_sat_array, Psat_array, label=r'$p_\mathrm{sat}$'+vol_latex[vol], lw=ls_ind, ls=":", color=vol_colors[vol][4]) + ax1.semilogy( T_sat_array, Psat_array, label=r'$p_\mathrm{sat}$'+vol_latex[vol], lw=ls_ind, ls=":", color=vol_colors[vol]) # Plot partial pressures - ax1.semilogy(atm.tmp, atm.p_vol[vol], color=vol_colors[vol][4], lw=ls_ind, ls="-", label=r'$p$'+vol_latex[vol],alpha=0.99) + ax1.semilogy(atm.tmp, atm.p_vol[vol], color=vol_colors[vol], lw=ls_ind, ls="-", label=r'$p$'+vol_latex[vol],alpha=0.99) # Sum up partial pressures p_partial_sum += atm.p_vol[vol] # Plot individual molar concentrations - ax2.semilogy(atm.x_cond[vol],atm.p, color=vol_colors[vol][4], lw=ls_ind, ls="--", label=vol_latex[vol]+" cond.") - ax2.semilogy(atm.x_gas[vol],atm.p, color=vol_colors[vol][4], lw=ls_ind, ls="-", label=vol_latex[vol]+" gas") + ax2.semilogy(atm.x_cond[vol],atm.p, color=vol_colors[vol], lw=ls_ind, ls="--", label=vol_latex[vol]+" cond.") + ax2.semilogy(atm.x_gas[vol],atm.p, color=vol_colors[vol], lw=ls_ind, ls="-", label=vol_latex[vol]+" gas") # # Plot sum of partial pressures as check ax1.semilogy(atm.tmp, p_partial_sum, color="green", lw=ls_dry, ls="--", label=r'$\sum p^\mathrm{i}$',alpha=0.99) diff --git a/src/janus/utils/data.py b/src/janus/utils/data.py index a1a4e32..6021a0e 100644 --- a/src/janus/utils/data.py +++ b/src/janus/utils/data.py @@ -1,17 +1,13 @@ import os from osfclient.api import OSF -#project ID of the stellar evolution tracks folder in the OSF -project_id = 'vehxg' - basic_list =[ "/Dayspring/256", "/Frostflow/256", "/Legacy", "/Mallard", "/Oak", - "/Reach", - "/stellar_spectra" + "/Reach" ] def download_folder(storage, folder_name, local_path): @@ -33,6 +29,31 @@ def download_folder(storage, folder_name, local_path): file.write_to(local_file) return + +def GetFWLData(): + fwl_data_dir = os.getenv('FWL_DATA') + if os.environ.get("FWL_DATA") == None: + raise Exception("The FWL_DATA environment variable where spectral data will be downloaded needs to be set up!") + return fwl_data_dir + +def DownloadStellarSpectra(): + #project ID of the stellar spectra on OSF + project_id = '8r2sw' + + # Link with OSF project repository + osf = OSF() + project = osf.project(project_id) + storage = project.storage('osfstorage') + + # Folder + data_dir = GetFWLData() + "/stellar_spectra" + if not os.path.exists(data_dir): + os.makedirs(data_dir) + + # Get all named spectra + download_folder(storage,"/Named",data_dir) + + def DownloadSpectralFiles(fname="",nband=256): '''' Download spectral files data @@ -44,13 +65,11 @@ def DownloadSpectralFiles(fname="",nband=256): (only relevant for Dayspring, Frostflow and Honeyside) ''' - #Check if data environment variable is set up - fwl_data_dir = os.getenv('FWL_DATA') - if os.environ.get("FWL_DATA") == None: - raise Exception("The FWL_DATA environment variable where spectral data will be downloaded needs to be set up!") + #project ID of the spectral files on OSF + project_id = 'vehxg' #Create spectral file data repository if not existing - data_dir = fwl_data_dir + "/spectral_files" + data_dir = GetFWLData() + "/spectral_files" if not os.path.exists(data_dir): os.makedirs(data_dir) diff --git a/src/janus/utils/phys.py b/src/janus/utils/phys.py index 1ba6b14..06776cf 100644 --- a/src/janus/utils/phys.py +++ b/src/janus/utils/phys.py @@ -76,7 +76,9 @@ "NO" : 0.0300061, "NO2" : 0.0460055, "HNO3": 0.06301284, - "OCS" : 0.060075 + "OCS" : 0.060075, + "S2" : 0.0641, # kg mol-1 + } #This class allows convenient access diff --git a/tests/test_instellation.py b/tests/test_instellation.py index 986efc4..39461c8 100755 --- a/tests/test_instellation.py +++ b/tests/test_instellation.py @@ -5,7 +5,7 @@ import numpy as np from janus.modules import MCPA_CBL -from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, ReadBandEdges, StellarSpectrum +from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, DownloadStellarSpectra, ReadBandEdges, StellarSpectrum import mors def test_instellation(): @@ -27,7 +27,7 @@ def test_instellation(): #Download required spectral files DownloadSpectralFiles("/Oak") - DownloadSpectralFiles("/stellar_spectra") + DownloadStellarSpectra() # Setup spectral file print("Inserting stellar spectrum") diff --git a/tests/test_runaway_greenhouse.py b/tests/test_runaway_greenhouse.py index 0906a63..812efb4 100755 --- a/tests/test_runaway_greenhouse.py +++ b/tests/test_runaway_greenhouse.py @@ -6,7 +6,7 @@ from importlib.resources import files from janus.modules import RadConvEqm -from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, ReadBandEdges, StellarSpectrum +from janus.utils import atmos, CleanOutputDir, DownloadSpectralFiles, DownloadStellarSpectra, ReadBandEdges, StellarSpectrum import mors def test_runaway_greenhouse(): @@ -28,7 +28,7 @@ def test_runaway_greenhouse(): #Download required spectral files DownloadSpectralFiles("/Oak") - DownloadSpectralFiles("/stellar_spectra") + DownloadStellarSpectra() # Setup spectral file print("Inserting stellar spectrum") From bafbffc2f732b4309b4a804f313a138bc33f1242 Mon Sep 17 00:00:00 2001 From: Harrison Nicholls Date: Thu, 18 Jul 2024 12:14:28 +0100 Subject: [PATCH 2/5] Fixes to download spectrum. SocRadConv works --- examples/SocRadConv.py | 10 +++++++++- examples/demo_instellation.py | 11 ++++++++++- examples/demo_runaway_greenhouse.py | 11 ++++++++++- src/janus/utils/__init__.py | 1 + src/janus/utils/data.py | 16 +++++++++++++--- tests/test_instellation.py | 10 +++++++++- tests/test_runaway_greenhouse.py | 10 +++++++++- 7 files changed, 61 insertions(+), 8 deletions(-) diff --git a/examples/SocRadConv.py b/examples/SocRadConv.py index 5b4cd40..84e4571 100755 --- a/examples/SocRadConv.py +++ b/examples/SocRadConv.py @@ -72,11 +72,19 @@ DownloadSpectralFiles("/Dayspring") DownloadStellarSpectra() + # Read spectrum + spec = mors.Spectrum() + spec.LoadTSV(os.environ.get('FWL_DATA')+"/stellar_spectra/Named/sun.txt") + + # Convert to SOCRATES format + socstar = os.path.join(dirs["output"], "socstar.txt") + StellarSpectrum.PrepareStellarSpectrum(spec.wl, spec.fl, socstar) + # Move/prepare spectral file print("Inserting stellar spectrum") StellarSpectrum.InsertStellarSpectrum( os.environ.get('FWL_DATA')+"/spectral_files/Dayspring/256/Dayspring.sf", - os.environ.get('FWL_DATA')+"/spectral_files/stellar_spectra/Sun_t4_4Ga_claire_12.txt", + socstar, dirs["output"] ) diff --git a/examples/demo_instellation.py b/examples/demo_instellation.py index 160b3d1..e498ec5 100755 --- a/examples/demo_instellation.py +++ b/examples/demo_instellation.py @@ -39,11 +39,20 @@ DownloadSpectralFiles("/Oak") DownloadStellarSpectra() + # Read spectrum + spec = mors.Spectrum() + spec.LoadTSV(os.environ.get('FWL_DATA')+"/stellar_spectra/Named/sun.txt") + + # Convert to SOCRATES format + socstar = os.path.join(dirs["output"], "socstar.txt") + StellarSpectrum.PrepareStellarSpectrum(spec.wl, spec.fl, socstar) + + # Setup spectral file print("Inserting stellar spectrum") StellarSpectrum.InsertStellarSpectrum( os.environ.get('FWL_DATA')+"/spectral_files/Oak/318/Oak.sf", - os.environ.get('FWL_DATA')+"/spectral_files/stellar_spectra/Sun_t4_4Ga_claire_12.txt", + socstar, dirs["output"] ) band_edges = ReadBandEdges(dirs["output"]+"star.sf") diff --git a/examples/demo_runaway_greenhouse.py b/examples/demo_runaway_greenhouse.py index 5b03165..47ab9d2 100755 --- a/examples/demo_runaway_greenhouse.py +++ b/examples/demo_runaway_greenhouse.py @@ -38,11 +38,20 @@ DownloadSpectralFiles("/Oak") DownloadStellarSpectra() + # Read spectrum + spec = mors.Spectrum() + spec.LoadTSV(os.environ.get('FWL_DATA')+"/stellar_spectra/Named/sun.txt") + + # Convert to SOCRATES format + socstar = os.path.join(dirs["output"], "socstar.txt") + StellarSpectrum.PrepareStellarSpectrum(spec.wl, spec.fl, socstar) + + # Setup spectral file print("Inserting stellar spectrum") StellarSpectrum.InsertStellarSpectrum( os.environ.get('FWL_DATA')+"/spectral_files/Oak/318/Oak.sf", - os.environ.get('FWL_DATA')+"/spectral_files/stellar_spectra/Sun_t4_4Ga_claire_12.txt", + socstar, dirs["output"] ) print(" ") diff --git a/src/janus/utils/__init__.py b/src/janus/utils/__init__.py index 7ffc2b9..c7e73f8 100644 --- a/src/janus/utils/__init__.py +++ b/src/janus/utils/__init__.py @@ -4,6 +4,7 @@ # Data download from .data import DownloadSpectralFiles +from .data import DownloadStellarSpectra # Socrates utility module from .socrates import CleanOutputDir diff --git a/src/janus/utils/data.py b/src/janus/utils/data.py index 6021a0e..4c77e3d 100644 --- a/src/janus/utils/data.py +++ b/src/janus/utils/data.py @@ -31,15 +31,22 @@ def download_folder(storage, folder_name, local_path): def GetFWLData(): + '''' + Get path to FWL data directory on the disk + ''' fwl_data_dir = os.getenv('FWL_DATA') if os.environ.get("FWL_DATA") == None: raise Exception("The FWL_DATA environment variable where spectral data will be downloaded needs to be set up!") - return fwl_data_dir + return os.path.abspath(fwl_data_dir) def DownloadStellarSpectra(): + '''' + Download stellar spectra + ''' + #project ID of the stellar spectra on OSF project_id = '8r2sw' - + # Link with OSF project repository osf = OSF() project = osf.project(project_id) @@ -51,7 +58,9 @@ def DownloadStellarSpectra(): os.makedirs(data_dir) # Get all named spectra - download_folder(storage,"/Named",data_dir) + if not os.path.exists(data_dir+"/Named"): + print("Downloading stellar spectra") + download_folder(storage,"/Named",data_dir) def DownloadSpectralFiles(fname="",nband=256): @@ -82,6 +91,7 @@ def DownloadSpectralFiles(fname="",nband=256): if not fname: for folder in basic_list: if not os.path.exists(data_dir+folder): + print("Downloading basic SOCRATES spectral files") download_folder(storage,folder,data_dir) elif fname in ["/Dayspring","/Frostflow","/Honeyside"]: folder = fname + "/" + str(nband) diff --git a/tests/test_instellation.py b/tests/test_instellation.py index 39461c8..593d09e 100755 --- a/tests/test_instellation.py +++ b/tests/test_instellation.py @@ -29,11 +29,19 @@ def test_instellation(): DownloadSpectralFiles("/Oak") DownloadStellarSpectra() + # Read spectrum + spec = mors.Spectrum() + spec.LoadTSV(os.environ.get('FWL_DATA')+"/stellar_spectra/Named/sun.txt") + + # Convert to SOCRATES format + socstar = os.path.join(dirs["output"], "socstar.txt") + StellarSpectrum.PrepareStellarSpectrum(spec.wl, spec.fl, socstar) + # Setup spectral file print("Inserting stellar spectrum") StellarSpectrum.InsertStellarSpectrum( os.environ.get('FWL_DATA')+"/spectral_files/Oak/318/Oak.sf", - os.environ.get('FWL_DATA')+"/spectral_files/stellar_spectra/Sun_t4_4Ga_claire_12.txt", + socstar, dirs["output"] ) band_edges = ReadBandEdges(dirs["output"]+"star.sf") diff --git a/tests/test_runaway_greenhouse.py b/tests/test_runaway_greenhouse.py index 812efb4..3b36a61 100755 --- a/tests/test_runaway_greenhouse.py +++ b/tests/test_runaway_greenhouse.py @@ -30,11 +30,19 @@ def test_runaway_greenhouse(): DownloadSpectralFiles("/Oak") DownloadStellarSpectra() + # Read spectrum + spec = mors.Spectrum() + spec.LoadTSV(os.environ.get('FWL_DATA')+"/stellar_spectra/Named/sun.txt") + + # Convert to SOCRATES format + socstar = os.path.join(dirs["output"], "socstar.txt") + StellarSpectrum.PrepareStellarSpectrum(spec.wl, spec.fl, socstar) + # Setup spectral file print("Inserting stellar spectrum") StellarSpectrum.InsertStellarSpectrum( os.environ.get('FWL_DATA')+"/spectral_files/Oak/318/Oak.sf", - os.environ.get('FWL_DATA')+"/spectral_files/stellar_spectra/Sun_t4_4Ga_claire_12.txt", + socstar, dirs["output"] ) band_edges = ReadBandEdges(dirs["output"]+"star.sf") From 5880a6fe13fe9f78e51f16f5208aafd0b33e706b Mon Sep 17 00:00:00 2001 From: Harrison Nicholls Date: Thu, 18 Jul 2024 13:38:24 +0100 Subject: [PATCH 3/5] Using interpolation to down-bin stellar spectra, since old method sometimes caused weird behaviour. Updated tests accordingly, with minor changes overall --- src/janus/data/tests/data_instellation.csv | 14 ++--- src/janus/data/tests/old.csv | 8 +++ src/janus/modules/plot_emission_spectrum.py | 12 +++- src/janus/utils/StellarSpectrum.py | 63 +++++---------------- src/janus/utils/socrates.py | 7 ++- tests/test_instellation.py | 24 +++++--- 6 files changed, 60 insertions(+), 68 deletions(-) create mode 100644 src/janus/data/tests/old.csv diff --git a/src/janus/data/tests/data_instellation.csv b/src/janus/data/tests/data_instellation.csv index 31fc236..5815bb1 100644 --- a/src/janus/data/tests/data_instellation.csv +++ b/src/janus/data/tests/data_instellation.csv @@ -1,8 +1,8 @@ # r [AU], S_0 [W m-2], OLR [W m-2], net [W m-2], ts [K], tr[K] -3.00000e-01,2.34296e+03,1.94531e+03,-7.41821e+01,3.00037e+03,4.19568e+02 -4.83333e-01,9.02638e+02,8.44987e+02,6.70043e+01,2.99966e+03,3.30552e+02 -6.66667e-01,4.74449e+02,5.02546e+02,9.36275e+01,2.99953e+03,2.81455e+02 -8.50000e-01,2.91856e+02,4.21740e+02,1.70194e+02,2.99915e+03,2.49261e+02 -1.03333e+00,1.97481e+02,4.17035e+02,2.46829e+02,2.99877e+03,2.26070e+02 -1.21667e+00,1.42450e+02,4.16624e+02,2.93848e+02,2.99853e+03,2.08342e+02 -1.40000e+00,1.07585e+02,4.16473e+02,3.23747e+02,2.99838e+03,1.94222e+02 +3.00000e-01,2.31540e+03,1.94626e+03,-3.69138e+02,3.00185e+03,4.19568e+02 +4.83333e-01,8.92021e+02,8.45389e+02,-4.66315e+01,3.00023e+03,3.30552e+02 +6.66667e-01,4.68869e+02,5.02606e+02,3.37376e+01,2.99983e+03,2.81454e+02 +8.50000e-01,2.88424e+02,4.21925e+02,1.33502e+02,2.99933e+03,2.49260e+02 +1.03333e+00,1.95159e+02,4.17176e+02,2.22017e+02,2.99889e+03,2.26070e+02 +1.21667e+00,1.40775e+02,4.16726e+02,2.75951e+02,2.99862e+03,2.08342e+02 +1.40000e+00,1.06319e+02,4.16550e+02,3.10230e+02,2.99845e+03,1.94222e+02 diff --git a/src/janus/data/tests/old.csv b/src/janus/data/tests/old.csv new file mode 100644 index 0000000..31fc236 --- /dev/null +++ b/src/janus/data/tests/old.csv @@ -0,0 +1,8 @@ +# r [AU], S_0 [W m-2], OLR [W m-2], net [W m-2], ts [K], tr[K] +3.00000e-01,2.34296e+03,1.94531e+03,-7.41821e+01,3.00037e+03,4.19568e+02 +4.83333e-01,9.02638e+02,8.44987e+02,6.70043e+01,2.99966e+03,3.30552e+02 +6.66667e-01,4.74449e+02,5.02546e+02,9.36275e+01,2.99953e+03,2.81455e+02 +8.50000e-01,2.91856e+02,4.21740e+02,1.70194e+02,2.99915e+03,2.49261e+02 +1.03333e+00,1.97481e+02,4.17035e+02,2.46829e+02,2.99877e+03,2.26070e+02 +1.21667e+00,1.42450e+02,4.16624e+02,2.93848e+02,2.99853e+03,2.08342e+02 +1.40000e+00,1.07585e+02,4.16473e+02,3.23747e+02,2.99838e+03,1.94222e+02 diff --git a/src/janus/modules/plot_emission_spectrum.py b/src/janus/modules/plot_emission_spectrum.py index e3ae616..a85be03 100644 --- a/src/janus/modules/plot_emission_spectrum.py +++ b/src/janus/modules/plot_emission_spectrum.py @@ -7,7 +7,8 @@ import janus.utils.phys as phys def plot_emission(atm:atmos, filename:str='output/toa_emission.pdf', - level_idx:int=0, planck_surface:bool=True, show_bands:bool=False): + level_idx:int=0, planck_surface:bool=True, + incoming_solar:bool=True, show_bands:bool=False): # Configure planck_samps = 1000 @@ -34,6 +35,13 @@ def plot_emission(atm:atmos, filename:str='output/toa_emission.pdf', # Make plot fig,ax = plt.subplots(1,1, figsize=(8,4)) + # Plot incoming stellar spectrum + if incoming_solar: + ysolar = atm.SW_spectral_flux_down.T[0] # W m-2 + ysolar /= atm.band_widths + ysolar *= 1e3 + ax.plot(x, ysolar, color="green", lw=lw, label="TOA stellar flux") + # Plot band edges if show_bands: for i in range(atm.nbands): @@ -70,7 +78,7 @@ def plot_emission(atm:atmos, filename:str='output/toa_emission.pdf', ax.set_xlim(max(xmin, np.amin(x)), min(xmax, np.amax(x))) ax.set_xscale("log") - ax.set_ylabel("Upward flux [erg s$^{-1}$ cm$^{-2}$ nm$^{-1}$]") + ax.set_ylabel("Flux density [erg s$^{-1}$ cm$^{-2}$ nm$^{-1}$]") ax.set_yscale("log") ax.legend(loc="lower center") diff --git a/src/janus/utils/StellarSpectrum.py b/src/janus/utils/StellarSpectrum.py index f33aa30..4173433 100644 --- a/src/janus/utils/StellarSpectrum.py +++ b/src/janus/utils/StellarSpectrum.py @@ -9,7 +9,7 @@ import numpy as np import shutil , os import subprocess -from scipy.stats import binned_statistic +from scipy.interpolate import PchipInterpolator def PrepareStellarSpectrum(wl, fl, star_file, nbins_max=95000): """Write a stellar spectrum. @@ -45,13 +45,9 @@ def PrepareStellarSpectrum(wl, fl, star_file, nbins_max=95000): raise Exception("Too many bins requested for stellar spectrum (maximum is %d)" % socrates_nbins_max) # Down-sample spectrum when necessary or requested - if (nbins_max > -1) or (len(wl) > socrates_nbins_max): + if len(wl) > socrates_nbins_max: print("Rebinning stellar spectrum") - # Parameters - max_retry = 10 # Number of times to down-sample before giving up - nbins_floor = 250 # Minimum number of bins - # Store old wl,fl arrays wl_orig = wl fl_orig = fl @@ -61,42 +57,13 @@ def PrepareStellarSpectrum(wl, fl, star_file, nbins_max=95000): nbins_max = min( int(socrates_nbins_max), nbins_max) # Must be fewer than 100k - # Perform rounds of down-sampling until nbins < nbins_max - # Rarely performs more than 2 rounds - bins_retry = True - downsample_factor = len(wl)/nbins_max - count_retry = 0 - while bins_retry: + # create interpolator on log-wavelength grid + itp = PchipInterpolator(np.log10(wl_orig), fl_orig) + wl = np.linspace(np.log10(wl_orig[0]), np.log10(wl_orig[-1]), nbins_max) - # Prevent loop from getting stuck - if count_retry > max_retry: - raise Exception("Giving up downsampling spectrum after %d rounds" % count_retry) - - # Ensure that the bins can be subdivided by truncating the array - downsample_factor = int(downsample_factor) - crop_amount = len(wl_orig) - int(len(wl_orig) % downsample_factor) - wl = wl_orig[:crop_amount] - - # Downsample - wl = wl.reshape(-1, downsample_factor).mean(axis=1) # https://saturncloud.io/blog/the-optimal-approach-to-downsampling-a-numpy-array/ - nbins = len(wl) - - if (nbins < nbins_floor): - raise Exception("Too few bins in stellar spectrum (%d bins)" % nbins) - else: - values, edges, _ = binned_statistic(wl_orig, fl_orig, statistic='mean', bins=wl) - wl = 0.5 * (edges[:-1] + edges[1:]) - fl = values - - if np.isnan(fl).any() or np.any(fl <= 0): - # Try again with fewer bins (occurs when there are empty bins) - downsample_factor *= 1.5 - else: - # Try again if still too many samples (unlikely case, but does happen sometimes) - bins_retry = bool( nbins > nbins_max ) - downsample_factor *= 1.2 - - count_retry += 1 + # interpolate new fluxes + fl = itp(wl) + wl = 10**wl # Convert units wl = np.array(wl) * 1.0e-9 # [nm] -> [m] @@ -156,14 +123,12 @@ def InsertStellarSpectrum(orig_file:str, star_file:str, output_folder:str): shutil.copyfile(orig_file, outp_file) shutil.copyfile(orig_filek, outp_filek) - # Copy star file to output folder - star_file_cpy = os.path.join(output_folder, "star.dat") - if os.path.exists(star_file_cpy): - os.remove(star_file_cpy) - shutil.copyfile(star_file, star_file_cpy) - - # Run prep_spec from SOCRATES to insert stellar spectrum - inputs = [outp_file,'a','6','n','T','100 4000','150','2','n',star_file_cpy,'y','-1','EOF'] + # Run prep_spec + inputs = [outp_file,'a', # append existing file + '6','n','T','100 4000','250', # tabulate thermal source function + '2','n',star_file,'y', # insert stellar spectrum + '-1','EOF' # save and exit + ] p = subprocess.run(['prep_spec'], stdout=subprocess.PIPE, input='\n'.join(inputs), encoding='ascii') if (p.returncode != 0): print("WARNING: prep_spec returned with code %d" % p.returncode) diff --git a/src/janus/utils/socrates.py b/src/janus/utils/socrates.py index 341edbf..2c3f8bb 100644 --- a/src/janus/utils/socrates.py +++ b/src/janus/utils/socrates.py @@ -52,7 +52,10 @@ def radCompSoc(atm, dirs, recalc, rscatter=False, # Only supported from SOCRATES version 2306 onwards (revision 1403) socrates_use_namelist = True - + # Change directory + olddir = os.getcwd() + os.chdir(dirs["output"]) + # Define path spectral files starspectral_file = dirs["output"]+"star.sf" runspectral_file = dirs["output"]+"runtime.sf" @@ -309,6 +312,8 @@ def radCompSoc(atm, dirs, recalc, rscatter=False, if atm.has_contfunc: ncfile11.close() + os.chdir(olddir) + return atm # Sting sorting not based on natsorted package diff --git a/tests/test_instellation.py b/tests/test_instellation.py index 593d09e..5898bb2 100755 --- a/tests/test_instellation.py +++ b/tests/test_instellation.py @@ -49,7 +49,7 @@ def test_instellation(): # Open config file cfg_file = dirs["janus"]+"data/tests/config_instellation.toml" with open(cfg_file, 'r') as f: - cfg = toml.load(f) + cfg = toml.load(f) # Star luminosity time = { "planet": cfg['planet']['time'], "star": cfg['star']['time']} @@ -73,17 +73,23 @@ def test_instellation(): r_arr = np.linspace(0.3, 1.4, 7) # orbital distance range [AU] for i in range(7): - print("Orbital separation = %.2f AU" % r_arr[i]) + print("Orbital separation = %.2f AU" % r_arr[i]) - atm.instellation = baraffe.BaraffeSolarConstant(time['star'], r_arr[i]) - atmos.setTropopauseTemperature(atm) + atm.instellation = baraffe.BaraffeSolarConstant(time['star'], r_arr[i]) + atmos.setTropopauseTemperature(atm) - atm = MCPA_CBL(dirs, atm, False, rscatter = True, T_surf_max=9.0e99, T_surf_guess = atm.trppT+100) + atm = MCPA_CBL(dirs, atm, False, rscatter = True, T_surf_max=9.0e99, T_surf_guess = atm.trppT+100) - out = [atm.SW_flux_down[0], atm.LW_flux_up[0], atm.net_flux[0], atm.ts, atm.trppT] - print(out) - print(ref[i][1:6]) - np.testing.assert_allclose(out, ref[i][1:6], rtol=1e-5, atol=0) + out = [atm.SW_flux_down[0], atm.LW_flux_up[0], atm.net_flux[0], atm.ts, atm.trppT] + + print_out = "%.5e,"%float(r_arr[i]) + print_out += ",".join(["%.5e"%o for o in out]) + print(print_out) + + print_out = "%.5e,"%float(r_arr[i]) + print_out += ",".join(["%.5e"%o for o in ref[i]]) + + np.testing.assert_allclose(out, ref[i][1:6], rtol=1e-5, atol=0) # Tidy CleanOutputDir(os.getcwd()) From b6ec61d77638707c230a723bfff3e1d55b322c29 Mon Sep 17 00:00:00 2001 From: Harrison Nicholls Date: Thu, 18 Jul 2024 14:43:16 +0100 Subject: [PATCH 4/5] Updated tests. Minor changes from master overall --- README.md | 7 ++-- src/janus/data/tests/data_instellation.csv | 14 +++---- .../data/tests/data_runaway_greenhouse.csv | 40 +++++++++---------- src/janus/data/tests/old.csv | 8 ---- tests/test_instellation.py | 6 +-- tests/test_runaway_greenhouse.py | 4 +- 6 files changed, 36 insertions(+), 43 deletions(-) delete mode 100644 src/janus/data/tests/old.csv diff --git a/README.md b/README.md index ac7065b..fd83b70 100644 --- a/README.md +++ b/README.md @@ -36,12 +36,13 @@ https://proteus-code.readthedocs.io * `git clone git@github.com:FormingWorlds/JANUS.git` * `cd JANUS` * `pip install -e .` -3. Download spectral files from the [OSF repository](https://osf.io/vehxg/) +3. Download data from the [OSF repository](https://osf.io/vehxg/) * Set the environment variable FWL_DATA to define where the spectral data files will be stored * `export FWL_DATA=...` - * Run the following commands within a python environment (or script) to download all basic spectral files - * `from janus.utils.data import DownloadSpectralFiles` + * Run the following commands within a python environment (or script) to download all basic data + * `from janus.utils.data import *` * `DownloadSpectralFiles()` + * `DownloadStellarSpectra()` * Alternatively, you can specify which spectral data you want to download, and optionally the number of bands * `DownloadSpectralFiles("/Frostflow", 4096)` diff --git a/src/janus/data/tests/data_instellation.csv b/src/janus/data/tests/data_instellation.csv index 5815bb1..c4aea42 100644 --- a/src/janus/data/tests/data_instellation.csv +++ b/src/janus/data/tests/data_instellation.csv @@ -1,8 +1,8 @@ # r [AU], S_0 [W m-2], OLR [W m-2], net [W m-2], ts [K], tr[K] -3.00000e-01,2.31540e+03,1.94626e+03,-3.69138e+02,3.00185e+03,4.19568e+02 -4.83333e-01,8.92021e+02,8.45389e+02,-4.66315e+01,3.00023e+03,3.30552e+02 -6.66667e-01,4.68869e+02,5.02606e+02,3.37376e+01,2.99983e+03,2.81454e+02 -8.50000e-01,2.88424e+02,4.21925e+02,1.33502e+02,2.99933e+03,2.49260e+02 -1.03333e+00,1.95159e+02,4.17176e+02,2.22017e+02,2.99889e+03,2.26070e+02 -1.21667e+00,1.40775e+02,4.16726e+02,2.75951e+02,2.99862e+03,2.08342e+02 -1.40000e+00,1.06319e+02,4.16550e+02,3.10230e+02,2.99845e+03,1.94222e+02 +3.00000e-01,2.34297e+03,1.94397e+03,-4.67185e+01,3.00023e+03,4.19568e+02 +4.83333e-01,9.02641e+02,8.43916e+02,7.70272e+01,2.99961e+03,3.30552e+02 +6.66667e-01,4.74451e+02,5.01731e+02,9.86425e+01,2.99951e+03,2.81454e+02 +8.50000e-01,2.91857e+02,4.21028e+02,1.73069e+02,2.99913e+03,2.49260e+02 +1.03333e+00,1.97482e+02,4.16319e+02,2.48540e+02,2.99876e+03,2.26070e+02 +1.21667e+00,1.42451e+02,4.15915e+02,2.94890e+02,2.99853e+03,2.08342e+02 +1.40000e+00,1.07585e+02,4.15768e+02,3.24365e+02,2.99838e+03,1.94222e+02 diff --git a/src/janus/data/tests/data_runaway_greenhouse.csv b/src/janus/data/tests/data_runaway_greenhouse.csv index b7c8618..724e9e2 100644 --- a/src/janus/data/tests/data_runaway_greenhouse.csv +++ b/src/janus/data/tests/data_runaway_greenhouse.csv @@ -1,21 +1,21 @@ # T_surf [K], OLR [W m-2] - 200.00000000, 90.72589 - 336.84210526, 277.55026 - 473.68421053, 277.47986 - 610.52631579, 277.27164 - 747.36842105, 277.35020 - 884.21052632, 277.18277 -1021.05263158, 277.37436 -1157.89473684, 277.13092 -1294.73684211, 277.14750 -1431.57894737, 277.30030 -1568.42105263, 277.64148 -1705.26315789, 279.04703 -1842.10526316, 285.55804 -1978.94736842, 308.2808 -2115.78947368, 375.7102 -2252.63157895, 548.4764 -2389.47368421, 954.2722 -2526.31578947, 1821.3682 -2663.15789474, 3517.7295 -2800.00000000, 6636.2110 + 200.00000000,9.07314e+01 + 336.84210526,2.77416e+02 + 473.68421053,2.77349e+02 + 610.52631579,2.77138e+02 + 747.36842105,2.77217e+02 + 884.21052632,2.77049e+02 +1021.05263158,2.77241e+02 +1157.89473684,2.76997e+02 +1294.73684211,2.77013e+02 +1431.57894737,2.77166e+02 +1568.42105263,2.77504e+02 +1705.26315789,2.78884e+02 +1842.10526316,2.85329e+02 +1978.94736842,3.07788e+02 +2115.78947368,3.74184e+02 +2252.63157895,5.47014e+02 +2389.47368421,9.50390e+02 +2526.31578947,1.80287e+03 +2663.15789474,3.49537e+03 +2800.00000000,6.58174e+03 diff --git a/src/janus/data/tests/old.csv b/src/janus/data/tests/old.csv deleted file mode 100644 index 31fc236..0000000 --- a/src/janus/data/tests/old.csv +++ /dev/null @@ -1,8 +0,0 @@ -# r [AU], S_0 [W m-2], OLR [W m-2], net [W m-2], ts [K], tr[K] -3.00000e-01,2.34296e+03,1.94531e+03,-7.41821e+01,3.00037e+03,4.19568e+02 -4.83333e-01,9.02638e+02,8.44987e+02,6.70043e+01,2.99966e+03,3.30552e+02 -6.66667e-01,4.74449e+02,5.02546e+02,9.36275e+01,2.99953e+03,2.81455e+02 -8.50000e-01,2.91856e+02,4.21740e+02,1.70194e+02,2.99915e+03,2.49261e+02 -1.03333e+00,1.97481e+02,4.17035e+02,2.46829e+02,2.99877e+03,2.26070e+02 -1.21667e+00,1.42450e+02,4.16624e+02,2.93848e+02,2.99853e+03,2.08342e+02 -1.40000e+00,1.07585e+02,4.16473e+02,3.23747e+02,2.99838e+03,1.94222e+02 diff --git a/tests/test_instellation.py b/tests/test_instellation.py index 5898bb2..e36c8f2 100755 --- a/tests/test_instellation.py +++ b/tests/test_instellation.py @@ -84,10 +84,10 @@ def test_instellation(): print_out = "%.5e,"%float(r_arr[i]) print_out += ",".join(["%.5e"%o for o in out]) - print(print_out) + print("Calculated:",print_out) - print_out = "%.5e,"%float(r_arr[i]) - print_out += ",".join(["%.5e"%o for o in ref[i]]) + # print_out = ",".join(["%.5e"%o for o in ref[i]]) + # print("Target:",print_out) np.testing.assert_allclose(out, ref[i][1:6], rtol=1e-5, atol=0) diff --git a/tests/test_runaway_greenhouse.py b/tests/test_runaway_greenhouse.py index 3b36a61..c81e124 100755 --- a/tests/test_runaway_greenhouse.py +++ b/tests/test_runaway_greenhouse.py @@ -83,8 +83,8 @@ def test_runaway_greenhouse(): _, atm_moist = RadConvEqm(dirs, time, atm, standalone=True, cp_dry=False, trppD=False, rscatter=False) - out = [atm_moist.LW_flux_up[0]] - print("Output %s; Reference %s" % (out, OLR_ref[i][1])) + out = atm_moist.LW_flux_up[0] + print("Output %.5e; Reference %.5e" % (out, OLR_ref[i][1])) np.testing.assert_allclose(out, OLR_ref[i][1], rtol=1e-5, atol=0) # Tidy From 57358ca552f81e0ac9cafb029e2a866b370caf92 Mon Sep 17 00:00:00 2001 From: Harrison Nicholls Date: Thu, 18 Jul 2024 15:20:32 +0100 Subject: [PATCH 5/5] Removed chdir socrates exec. Already implemented in PROTEUS --- src/janus/utils/socrates.py | 8 +------- 1 file changed, 1 insertion(+), 7 deletions(-) diff --git a/src/janus/utils/socrates.py b/src/janus/utils/socrates.py index 2c3f8bb..9dcdeb0 100644 --- a/src/janus/utils/socrates.py +++ b/src/janus/utils/socrates.py @@ -19,7 +19,7 @@ def radCompSoc(atm, dirs, recalc, rscatter=False, - rewrite_cfg=True, rewrite_tmp=True, rewrite_gas=True): + rewrite_cfg=True, rewrite_tmp=True, rewrite_gas=False): """Runs SOCRATES to calculate fluxes and heating rates Parameters @@ -52,10 +52,6 @@ def radCompSoc(atm, dirs, recalc, rscatter=False, # Only supported from SOCRATES version 2306 onwards (revision 1403) socrates_use_namelist = True - # Change directory - olddir = os.getcwd() - os.chdir(dirs["output"]) - # Define path spectral files starspectral_file = dirs["output"]+"star.sf" runspectral_file = dirs["output"]+"runtime.sf" @@ -312,8 +308,6 @@ def radCompSoc(atm, dirs, recalc, rscatter=False, if atm.has_contfunc: ncfile11.close() - os.chdir(olddir) - return atm # Sting sorting not based on natsorted package