Fixing ice radius/diameter bug and a whole bunch of other changes

allsky/allsky_init.py

@@ -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
 
@@ -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
@@ -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()

include/Cloud_optics.h

@@ -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:
@@ -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;
@@ -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:
@@ -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;

include_rt/Cloud_optics_rt.h

@@ -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:
@@ -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;

include_test/Radiation_solver.h

@@ -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,
@@ -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,
@@ -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,
@@ -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,

include_test/Radiation_solver_bw.h

@@ -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,
@@ -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,
@@ -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,
@@ -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,
@@ -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,
@@ -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

include_test/Radiation_solver_rt.h

@@ -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,
@@ -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,
@@ -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,

rfmip/rfmip_init.py

@@ -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()
 
@@ -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)
@@ -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)
@@ -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()