Skip to content

Commit

Permalink
Merge pull request #45 from MennoVeerman/main
Browse files Browse the repository at this point in the history
Fixing ice radius/diameter bug and a whole bunch of other changes
  • Loading branch information
Chiil authored Dec 19, 2024
2 parents db5dc0b + c1bb890 commit a0f96ac
Show file tree
Hide file tree
Showing 20 changed files with 334 additions and 309 deletions.
16 changes: 8 additions & 8 deletions allsky/allsky_init.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,7 @@ def calc_p_q_T(z):

i_above_zt = np.where(z > z_trop)
q[i_above_zt] = q_t

T_0 = 300.
Tv_0 = (1. + 0.608*q_0)*T_0

Expand All @@ -51,16 +51,16 @@ def calc_p_q_T(z):
Tv = T * (1. + 0.608*q)

# T = Tv / (1. + 0.608*q)

g = 9.79764
Rd = 287.04
p0 = 101480.

p = p0 * (Tv / Tv_0)**(g/(Rd*gamma))

p_tmp = p0 * (Tv/Tv_0)**(g/(Rd*gamma)) \
* np.exp( -( (g*(z-z_trop)) / (Rd*Tv) ) )

p[i_above_zt] = p_tmp[i_above_zt]

return p, q, T
Expand Down Expand Up @@ -159,20 +159,20 @@ def calc_p_q_T(z):
nc_lwp = nc_file.createVariable('lwp', float_type, ('lay', 'y', 'x'))
nc_iwp = nc_file.createVariable('iwp', float_type, ('lay', 'y', 'x'))
nc_rel = nc_file.createVariable('rel', float_type, ('lay', 'y', 'x'))
nc_rei = nc_file.createVariable('rei', float_type, ('lay', 'y', 'x'))
nc_dei = nc_file.createVariable('dei', float_type, ('lay', 'y', 'x'))

min_rel, max_rel = 2.5, 21.5
min_rei, max_rei = 10., 180.
min_dei, max_dei = 10., 180.

rel_val = 0.5*(min_rel + max_rel)
rei_val = 0.5*(min_rei + max_rei)
dei_val = 0.5*(min_dei + max_dei)

cloud_flag = (np.arange(1, n_col_x+1)%3 > 0)
cloud_mask = np.where((nc_p_lay[:,:,:] > 1.e4) & (nc_p_lay[:,:,:] < 9.e4) & cloud_flag[None, None,:], True, False)

nc_lwp[:,:,:] = np.where(cloud_mask & (nc_T_lay[:,:,:] > 263.), 10., 0.)
nc_iwp[:,:,:] = np.where(cloud_mask & (nc_T_lay[:,:,:] < 273.), 10., 0.)
nc_rel[:,:,:] = np.where(nc_lwp[:,:,:] > 0., rel_val, 0.)
nc_rei[:,:,:] = np.where(nc_iwp[:,:,:] > 0., rei_val, 0.)
nc_dei[:,:,:] = np.where(nc_iwp[:,:,:] > 0., dei_val, 0.)

nc_file.close()
20 changes: 10 additions & 10 deletions include/Cloud_optics.h
Original file line number Diff line number Diff line change
Expand Up @@ -41,18 +41,18 @@ class Cloud_optics : public Optical_props
Cloud_optics(
const Array<Float,2>& band_lims_wvn,
const Float radliq_lwr, const Float radliq_upr, const Float radliq_fac,
const Float radice_lwr, const Float radice_upr, const Float radice_fac,
const Float diamice_lwr, const Float diamice_upr, const Float diamice_fac,
const Array<Float,2>& lut_extliq, const Array<Float,2>& lut_ssaliq, const Array<Float,2>& lut_asyliq,
const Array<Float,3>& lut_extice, const Array<Float,3>& lut_ssaice, const Array<Float,3>& lut_asyice);

void cloud_optics(
const Array<Float,2>& clwp, const Array<Float,2>& ciwp,
const Array<Float,2>& reliq, const Array<Float,2>& reice,
const Array<Float,2>& reliq, const Array<Float,2>& deice,
Optical_props_1scl& optical_props);

void cloud_optics(
const Array<Float,2>& clwp, const Array<Float,2>& ciwp,
const Array<Float,2>& reliq, const Array<Float,2>& reice,
const Array<Float,2>& reliq, const Array<Float,2>& deice,
Optical_props_2str& optical_props);

private:
Expand All @@ -64,8 +64,8 @@ class Cloud_optics : public Optical_props
// Lookup table constants.
Float radliq_lwr;
Float radliq_upr;
Float radice_lwr;
Float radice_upr;
Float diamice_lwr;
Float diamice_upr;

// Lookup table coefficients.
Array<Float,2> lut_extliq;
Expand All @@ -84,18 +84,18 @@ class Cloud_optics_gpu : public Optical_props_gpu
Cloud_optics_gpu(
const Array<Float,2>& band_lims_wvn,
const Float radliq_lwr, const Float radliq_upr, const Float radliq_fac,
const Float radice_lwr, const Float radice_upr, const Float radice_fac,
const Float diamice_lwr, const Float diamice_upr, const Float diamice_fac,
const Array<Float,2>& lut_extliq, const Array<Float,2>& lut_ssaliq, const Array<Float,2>& lut_asyliq,
const Array<Float,3>& lut_extice, const Array<Float,3>& lut_ssaice, const Array<Float,3>& lut_asyice);

void cloud_optics(
const Array_gpu<Float,2>& clwp, const Array_gpu<Float,2>& ciwp,
const Array_gpu<Float,2>& reliq, const Array_gpu<Float,2>& reice,
const Array_gpu<Float,2>& reliq, const Array_gpu<Float,2>& deice,
Optical_props_1scl_gpu& optical_props);

void cloud_optics(
const Array_gpu<Float,2>& clwp, const Array_gpu<Float,2>& ciwp,
const Array_gpu<Float,2>& reliq, const Array_gpu<Float,2>& reice,
const Array_gpu<Float,2>& reliq, const Array_gpu<Float,2>& deice,
Optical_props_2str_gpu& optical_props);

private:
Expand All @@ -107,8 +107,8 @@ class Cloud_optics_gpu : public Optical_props_gpu
// Lookup table constants.
Float radliq_lwr;
Float radliq_upr;
Float radice_lwr;
Float radice_upr;
Float diamice_lwr;
Float diamice_upr;

// Lookup table coefficients.
Array<Float,2> lut_extliq;
Expand Down
12 changes: 6 additions & 6 deletions include_rt/Cloud_optics_rt.h
Original file line number Diff line number Diff line change
Expand Up @@ -32,27 +32,27 @@
// Forward declarations.
class Optical_props_rt;

#ifdef USECUDA
#ifdef USECUDA
class Cloud_optics_rt : public Optical_props_rt
{
public:
Cloud_optics_rt(
const Array<Float,2>& band_lims_wvn,
const Float radliq_lwr, const Float radliq_upr, const Float radliq_fac,
const Float radice_lwr, const Float radice_upr, const Float radice_fac,
const Float diamice_lwr, const Float diamice_upr, const Float diamice_fac,
const Array<Float,2>& lut_extliq, const Array<Float,2>& lut_ssaliq, const Array<Float,2>& lut_asyliq,
const Array<Float,3>& lut_extice, const Array<Float,3>& lut_ssaice, const Array<Float,3>& lut_asyice);

void cloud_optics(
const int ibnd,
const Array_gpu<Float,2>& clwp, const Array_gpu<Float,2>& ciwp,
const Array_gpu<Float,2>& reliq, const Array_gpu<Float,2>& reice,
const Array_gpu<Float,2>& reliq, const Array_gpu<Float,2>& deice,
Optical_props_1scl_rt& optical_props);

void cloud_optics(
const int ibnd,
const Array_gpu<Float,2>& clwp, const Array_gpu<Float,2>& ciwp,
const Array_gpu<Float,2>& reliq, const Array_gpu<Float,2>& reice,
const Array_gpu<Float,2>& reliq, const Array_gpu<Float,2>& deice,
Optical_props_2str_rt& optical_props);

private:
Expand All @@ -64,8 +64,8 @@ class Cloud_optics_rt : public Optical_props_rt
// Lookup table constants.
Float radliq_lwr;
Float radliq_upr;
Float radice_lwr;
Float radice_upr;
Float diamice_lwr;
Float diamice_upr;

// Lookup table coefficients.
Array<Float,2> lut_extliq;
Expand Down
8 changes: 4 additions & 4 deletions include_test/Radiation_solver.h
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,7 @@ class Radiation_solver_longwave
const Array<Float,2>& col_dry,
const Array<Float,1>& t_sfc, const Array<Float,2>& emis_sfc,
const Array<Float,2>& lwp, const Array<Float,2>& iwp,
const Array<Float,2>& rel, const Array<Float,2>& rei,
const Array<Float,2>& rel, const Array<Float,2>& dei,
Array<Float,3>& tau, Array<Float,3>& lay_source,
Array<Float,3>& lev_source, Array<Float,2>& sfc_source,
Array<Float,2>& lw_flux_up, Array<Float,2>& lw_flux_dn, Array<Float,2>& lw_flux_net,
Expand All @@ -81,7 +81,7 @@ class Radiation_solver_longwave
const Array_gpu<Float,2>& col_dry,
const Array_gpu<Float,1>& t_sfc, const Array_gpu<Float,2>& emis_sfc,
const Array_gpu<Float,2>& lwp, const Array_gpu<Float,2>& iwp,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& rei,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& dei,
Array_gpu<Float,3>& tau, Array_gpu<Float,3>& lay_source,
Array_gpu<Float,3>& lev_source, Array_gpu<Float,2>& sfc_source,
Array_gpu<Float,2>& lw_flux_up, Array_gpu<Float,2>& lw_flux_dn, Array_gpu<Float,2>& lw_flux_net,
Expand Down Expand Up @@ -151,7 +151,7 @@ class Radiation_solver_shortwave
const Array<Float,2>& sfc_alb_dir, const Array<Float,2>& sfc_alb_dif,
const Array<Float,1>& tsi_scaling, const Array<Float,1>& mu0,
const Array<Float,2>& lwp, const Array<Float,2>& iwp,
const Array<Float,2>& rel, const Array<Float,2>& rei,
const Array<Float,2>& rel, const Array<Float,2>& dei,
const Array<Float,2>& rh,
const Aerosol_concs& aerosol_concs,
Array<Float,3>& tau, Array<Float,3>& ssa, Array<Float,3>& g,
Expand Down Expand Up @@ -188,7 +188,7 @@ class Radiation_solver_shortwave
const Array_gpu<Float,2>& sfc_alb_dir, const Array_gpu<Float,2>& sfc_alb_dif,
const Array_gpu<Float,1>& tsi_scaling, const Array_gpu<Float,1>& mu0,
const Array_gpu<Float,2>& lwp, const Array_gpu<Float,2>& iwp,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& rei,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& dei,
const Array_gpu<Float,2>& rh,
const Aerosol_concs_gpu& aerosol_concs,
Array_gpu<Float,3>& tau, Array_gpu<Float,3>& ssa, Array_gpu<Float,3>& g,
Expand Down
14 changes: 7 additions & 7 deletions include_test/Radiation_solver_bw.h
Original file line number Diff line number Diff line change
Expand Up @@ -58,7 +58,7 @@ class Radiation_solver_longwave
const Array<Float,2>& col_dry,
const Array<Float,1>& t_sfc, const Array<Float,2>& emis_sfc,
const Array<Float,2>& lwp, const Array<Float,2>& iwp,
const Array<Float,2>& rel, const Array<Float,2>& rei,
const Array<Float,2>& rel, const Array<Float,2>& dei,
Array<Float,3>& tau, Array<Float,3>& lay_source,
Array<Float,3>& lev_source_inc, Array<Float,3>& lev_source_dec, Array<Float,2>& sfc_source,
Array<Float,2>& lw_flux_up, Array<Float,2>& lw_flux_dn, Array<Float,2>& lw_flux_net,
Expand All @@ -76,7 +76,7 @@ class Radiation_solver_longwave
Array_gpu<Float,2>& col_dry,
const Array_gpu<Float,1>& t_sfc, const Array_gpu<Float,2>& emis_sfc,
const Array_gpu<Float,2>& lwp, const Array_gpu<Float,2>& iwp,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& rei,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& dei,
Array_gpu<Float,3>& tau, Array_gpu<Float,3>& lay_source,
Array_gpu<Float,3>& lev_source_inc, Array_gpu<Float,3>& lev_source_dec, Array_gpu<Float,2>& sfc_source,
Array_gpu<Float,2>& lw_flux_up, Array_gpu<Float,2>& lw_flux_dn, Array_gpu<Float,2>& lw_flux_net,
Expand Down Expand Up @@ -136,7 +136,7 @@ class Radiation_solver_shortwave
const Array<Float,2>& sfc_alb,
const Array<Float,1>& tsi_scaling, const Array<Float,1>& mu0,
const Array<Float,2>& lwp, const Array<Float,2>& iwp,
const Array<Float,2>& rel, const Array<Float,2>& rei,
const Array<Float,2>& rel, const Array<Float,2>& dei,
Array<Float,3>& tau, Array<Float,3>& ssa, Array<Float,3>& g,
Array<Float,2>& toa_src,
Array<Float,2>& sw_flux_up, Array<Float,2>& sw_flux_dn,
Expand Down Expand Up @@ -174,7 +174,7 @@ class Radiation_solver_shortwave
const Array_gpu<Float,1>& tsi_scaling,
const Array_gpu<Float,1>& mu0, const Array_gpu<Float,1>& azi,
const Array_gpu<Float,2>& lwp, const Array_gpu<Float,2>& iwp,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& rei,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& dei,
const Array_gpu<Float,1>& land_use_map,
const Array_gpu<Float,2>& rh,
const Aerosol_concs_gpu& aerosol_concs,
Expand Down Expand Up @@ -209,7 +209,7 @@ class Radiation_solver_shortwave
const Array_gpu<Float,1>& tsi_scaling,
const Array_gpu<Float,1>& mu0, const Array_gpu<Float,1>& azi,
const Array_gpu<Float,2>& lwp, const Array_gpu<Float,2>& iwp,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& rei,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& dei,
const Array_gpu<Float,1>& land_use_map,
const Array_gpu<Float,2>& rh,
const Aerosol_concs_gpu& aerosol_concs,
Expand Down Expand Up @@ -250,12 +250,12 @@ class Radiation_solver_shortwave
Array_gpu<Float,3> mie_cdfs_bb;
Array_gpu<Float,4> mie_angs_bb;
Array_gpu<Float,4> mie_phase_bb;

Array_gpu<Float,1> mie_phase_angs;
Array_gpu<Float,3> mie_cdfs;
Array_gpu<Float,4> mie_angs;
Array_gpu<Float,4> mie_phase;

#endif
};
#endif
10 changes: 5 additions & 5 deletions include_test/Radiation_solver_rt.h
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,7 @@ class Radiation_solver_longwave
Array_gpu<Float,2>& col_dry,
const Array_gpu<Float,1>& t_sfc, const Array_gpu<Float,2>& emis_sfc,
const Array_gpu<Float,2>& lwp, const Array_gpu<Float,2>& iwp,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& rei,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& dei,
Array_gpu<Float,2>& tau, Array_gpu<Float,2>& lay_source,
Array_gpu<Float,2>& lev_source_inc, Array_gpu<Float,2>& lev_source_dec, Array_gpu<Float,1>& sfc_source,
Array_gpu<Float,2>& lw_flux_up, Array_gpu<Float,2>& lw_flux_dn, Array_gpu<Float,2>& lw_flux_net,
Expand Down Expand Up @@ -99,7 +99,7 @@ class Radiation_solver_shortwave
#ifdef __CUDACC__
void solve_gpu(
const bool switch_fluxes,
const bool switch_disable_2s,
const bool switch_twostream,
const bool switch_raytracing,
const bool switch_independent_column,
const bool switch_cloud_optics,
Expand All @@ -121,11 +121,11 @@ class Radiation_solver_shortwave
const Array_gpu<Float,1>& tsi_scaling,
const Array_gpu<Float,1>& mu0, const Array_gpu<Float,1>& azi,
const Array_gpu<Float,2>& lwp, const Array_gpu<Float,2>& iwp,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& rei,
const Array_gpu<Float,2>& rel, const Array_gpu<Float,2>& dei,
const Array_gpu<Float,2>& rh,
const Aerosol_concs_gpu& aerosol_concs,
Array_gpu<Float,2>& tot_tau_out, Array_gpu<Float,2>& tot_ssa_out,
Array_gpu<Float,2>& cld_tau_out, Array_gpu<Float,2>& cld_ssa_out, Array_gpu<Float,2>& cld_asy_out,
Array_gpu<Float,2>& tot_tau_out, Array_gpu<Float,2>& tot_ssa_out,
Array_gpu<Float,2>& cld_tau_out, Array_gpu<Float,2>& cld_ssa_out, Array_gpu<Float,2>& cld_asy_out,
Array_gpu<Float,2>& aer_tau_out, Array_gpu<Float,2>& aer_ssa_out, Array_gpu<Float,2>& aer_asy_out,
Array_gpu<Float,2>& sw_flux_up, Array_gpu<Float,2>& sw_flux_dn,
Array_gpu<Float,2>& sw_flux_dn_dir, Array_gpu<Float,2>& sw_flux_net,
Expand Down
28 changes: 14 additions & 14 deletions rfmip/rfmip_init.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,20 +12,20 @@
# Save all the input data to NetCDF
nc_file = nc.Dataset('rte_rrtmgp_input_expt_{:02d}.nc'.format(expt), mode='w', datamodel='NETCDF4', clobber=True)
nc_file_rfmip = nc.Dataset('multiple_input4MIPs_radiation_RFMIP_UColorado-RFMIP-1-2_none.nc', mode='r', datamodel='NETCDF4', clobber=False)

# Create a group for the radiation and set up the values.
nc_file.createDimension('lay', nc_file_rfmip.dimensions['layer'].size)
nc_file.createDimension('lev', nc_file_rfmip.dimensions['level'].size)
nc_file.createDimension('y', 1)
nc_file.createDimension('x', nc_file_rfmip.dimensions['site'].size)
nc_file.createDimension('band_lw', band_lw)
nc_file.createDimension('band_sw', band_sw)

nc_pres_layer = nc_file.createVariable('p_lay', float_type, ('lay', 'y', 'x'))
nc_pres_level = nc_file.createVariable('p_lev', float_type, ('lev', 'y', 'x'))
nc_temp_layer = nc_file.createVariable('t_lay', float_type, ('lay', 'y', 'x'))
nc_temp_level = nc_file.createVariable('t_lev', float_type, ('lev', 'y', 'x'))

nc_pres_layer[:,:,:] = nc_file_rfmip.variables['pres_layer'][:,:].transpose()
nc_pres_level[:,:,:] = nc_file_rfmip.variables['pres_level'][:,:].transpose()

Expand All @@ -35,25 +35,25 @@

nc_temp_layer[:,:,:] = (nc_file_rfmip.variables['temp_layer'][expt,:,:]).transpose()
nc_temp_level[:,:,:] = (nc_file_rfmip.variables['temp_level'][expt,:,:]).transpose()

nc_surface_emissivity = nc_file.createVariable('emis_sfc', float_type, ('y', 'x', 'band_lw'))
nc_surface_temperature = nc_file.createVariable('t_sfc', float_type, ('y', 'x'))

nc_surface_emissivity[:,:] = np.tile( (nc_file_rfmip.variables['surface_emissivity'][:]) [:,None], (1, band_lw) )
nc_surface_temperature[:,:] = nc_file_rfmip.variables['surface_temperature'][expt,:]

nc_surface_albedo_dir = nc_file.createVariable('sfc_alb_dir', float_type, ('y', 'x', 'band_sw'))
nc_surface_albedo_dif = nc_file.createVariable('sfc_alb_dif', float_type, ('y', 'x', 'band_sw'))

nc_surface_albedo_dir[:,:,:] = np.tile( (nc_file_rfmip.variables['surface_albedo'][:]) [:,None], (1, band_sw) )
nc_surface_albedo_dif[:,:,:] = np.tile( (nc_file_rfmip.variables['surface_albedo'][:]) [:,None], (1, band_sw) )

nc_mu0 = nc_file.createVariable('mu0', float_type, ('y', 'x'))
nc_mu0[:,:] = np.maximum(0., np.cos( np.deg2rad(nc_file_rfmip.variables['solar_zenith_angle'][:]) ) )

nc_tsi = nc_file.createVariable('tsi', float_type, ('y', 'x'))
nc_tsi[:,:] = nc_file_rfmip.variables['total_solar_irradiance'][:]

nc_h2o = nc_file.createVariable('vmr_h2o' , float_type, ('lay', 'y', 'x'))
nc_o3 = nc_file.createVariable('vmr_o3' , float_type, ('lay', 'y', 'x'))
nc_co2 = nc_file.createVariable('vmr_co2' , float_type)
Expand All @@ -77,7 +77,7 @@
# Profiles
nc_h2o[:,:,:] = nc_file_rfmip.variables['water_vapor'][expt,:,:].transpose() * float(nc_file_rfmip.variables['water_vapor'].units)
nc_o3 [:,:,:] = nc_file_rfmip.variables['ozone'][expt,:,:].transpose() * float(nc_file_rfmip.variables['ozone'].units)

# Constants
nc_co2 [:] = nc_file_rfmip.variables['carbon_dioxide_GM'][expt] * float(nc_file_rfmip.variables['carbon_dioxide_GM'].units)
nc_n2o [:] = nc_file_rfmip.variables['nitrous_oxide_GM'][expt] * float(nc_file_rfmip.variables['nitrous_oxide_GM'].units)
Expand All @@ -100,13 +100,13 @@
# CvH: To be removed if settings can be set.
nc_lwp = nc_file.createVariable('lwp', float_type, ('lay', 'y', 'x'))
nc_iwp = nc_file.createVariable('iwp', float_type, ('lay', 'y', 'x'))
nc_rel = nc_file.createVariable('rei', float_type, ('lay', 'y', 'x'))
nc_rei = nc_file.createVariable('rel', float_type, ('lay', 'y', 'x'))
nc_rel = nc_file.createVariable('dei', float_type, ('lay', 'y', 'x'))
nc_dei = nc_file.createVariable('rel', float_type, ('lay', 'y', 'x'))

nc_lwp[:,:,:] = 0.
nc_iwp[:,:,:] = 0.
nc_rel[:,:,:] = 0.
nc_rei[:,:,:] = 0.
nc_dei[:,:,:] = 0.
# CvH end.

nc_file_rfmip.close()
Expand Down
Loading

0 comments on commit a0f96ac

Please sign in to comment.