diff --git a/Source/Radiation/Albedo.H b/Source/Radiation/Albedo.H
new file mode 100644
index 000000000..bbb971fc8
--- /dev/null
+++ b/Source/Radiation/Albedo.H
@@ -0,0 +1,22 @@
+/*------------------------------------------------------------
+   Computes surface albedos
+   For more details , see Briegleb, Bruce P., 1992: Delta-Eddington
+   Approximation for Solar Radiation in the NCAR Community Climate Model,
+   Journal of Geophysical Research, Vol 97, D7, pp7603-7612).
+
+   NOTE: this is the simplest formula for albedo calculation, we have to
+         update to the latest implementation from microphysics output.
+*/
+#ifndef ERF_ALBEDO_H_
+#define ERF_ALBEDO_H_
+void set_albedo(const real1d& coszrs, real2d& albedo_dir, real2d& albedo_dif)
+{
+  // Albedos for land type I (Briegleb)
+  auto nswbands = albedo_dir.extent(0);
+  auto ncol     = albedo_dif.extent(1);
+  yakl::c::parallel_for(yakl::c::Bounds<2>(nswbands,ncol), YAKL_LAMBDA (int ibnd, int icol) {
+    albedo_dir(ibnd, icol) = 1.4 * 0.24 / ( 1. + 0.8 * coszrs(icol));
+    albedo_dif(ibnd, icol) = 1.2 * 0.24;
+  });
+}
+#endif
diff --git a/Source/Radiation/Init_rrtmgp.cpp b/Source/Radiation/Init_rrtmgp.cpp
index 12da46121..f00045cf4 100644
--- a/Source/Radiation/Init_rrtmgp.cpp
+++ b/Source/Radiation/Init_rrtmgp.cpp
@@ -10,7 +10,9 @@
 // Rrtmgp
 #include "Rrtmgp.H"
 
-void Rrtmgp::initialize()
+void Rrtmgp::initialize(int num_gas, std::vector<std::string> active_gas_names,
+                        const char* rrtmgp_coefficients_file_sw,
+                        const char* rrtmgp_coefficients_file_lw)
 {
     // First, make sure yakl has been initialized
     if (!yakl::isInitialized()) {
@@ -26,6 +28,13 @@ void Rrtmgp::initialize()
     // impossible from this initialization routine because I do not think the
     // rad_cnst objects are setup yet.
     // the other tasks!
+    ngas      = num_gas;
+    for (auto i = 0; i < ngas; ++i)
+       gas_names[i] = active_gas_names[i].c_str();
+
+    coefficients_file_sw = rrtmgp_coefficients_file_sw;
+    coefficients_file_lw = rrtmgp_coefficients_file_lw;
+
     auto active_gases = string1d("active_gases", ngas);
     for (int igas=0; igas<ngas; igas++) {
         active_gases(igas+1) = gas_names[igas];
diff --git a/Source/Radiation/Parameterizations.H b/Source/Radiation/Parameterizations.H
index 2131e04bf..c7c914fc5 100644
--- a/Source/Radiation/Parameterizations.H
+++ b/Source/Radiation/Parameterizations.H
@@ -39,67 +39,4 @@
     return RadConstants::retab[index]*(1.-fraction)+RadConstants::retab[index+1]*fraction;
  }
 
-  // Computes surface albedos over ocean and the surface
-  // Two spectral surface albedos for direct (dir) and diffuse (dif)
-  // incident radiation are calculated. The spectral intervals are:
-  //   s (shortwave)  = 0.2-0.7 micro-meters
-  //   l (longwave)   = 0.7-5.0 micro-meters
-  //
-  // Uses knowledge of surface type to specify albedo, as follows:
-  //
-  // Ocean           Uses solar zenith angle to compute albedo for direct
-  //                 radiation; diffuse radiation values constant; albedo
-  //                 independent of spectral interval and other physical
-  //                 factors such as ocean surface wind speed.
-  //
-  // For more details , see Briegleb, Bruce P., 1992: Delta-Eddington
-  // Approximation for Solar Radiation in the NCAR Community Climate Model,
-  // Journal of Geophysical Research, Vol 97, D7, pp7603-7612).
- AMREX_FORCE_INLINE
- void albedo(bool ocean, real1d& coszrs, real1d& tsurface, real1d& asdir,
-              real1d& aldir, real1d& asdif, real1d& aldif)
-  {
-    using yakl::intrinsics::size;
-    const real adif = 0.06;
-    int ncol = size(tsurface, 1);
-    if (ocean) {
-      // Ice-free ocean albedos function of solar zenith angle only, and
-      // independent of spectral interval:
-      yakl::c::parallel_for(ncol, YAKL_LAMBDA (int icol) {
-        if (coszrs(icol) <= 0) {
-          aldir(icol) = 0;
-          asdir(icol) = 0;
-          aldif(icol) = 0;
-          asdif(icol) = 0;
-        } else if (tsurface(icol) > 271.) {
-          aldir(icol) = ( .026 / (std::pow(coszrs(icol), 1.7) + .065)) +
-             (.15*(coszrs(icol) - 0.10) * (coszrs(icol) - 0.50) * (coszrs(icol) - 1.00) );
-          asdir(icol) = aldir(icol);
-          aldif(icol) = adif;
-          asdif(icol) = adif;
-        } else {
-         // albedo of sea-ice/snow surface
-           aldir(icol) = 0.45;
-           asdir(icol) = 0.75;
-           aldif(icol) = 0.45;
-           asdif(icol) = 0.75;
-       }
-      });
-    } else  {// land
-       yakl::c::parallel_for(ncol, YAKL_LAMBDA (int icol) {
-         if (coszrs(icol) <= 0.) {
-           aldir(icol) = 0.;
-           asdir(icol) = 0.;
-           aldif(icol) = 0.;
-           asdif(icol) = 0.;
-         } else {
-           // Albedos for land type I (Briegleb)
-           asdir(icol) = 1.4 * 0.06 / ( 1. + 0.8 * coszrs(icol));
-           asdif(icol) = 1.2 * 0.06;
-           aldir(icol) = 1.4 * 0.24 / ( 1. + 0.8 * coszrs(icol));
-           aldif(icol) = 1.2 * 0.24;
-        }
-     });
-   }
- }
 #endif
diff --git a/Source/Radiation/Rrtmgp.H b/Source/Radiation/Rrtmgp.H
index 8d0790f5b..9cca6a6f4 100644
--- a/Source/Radiation/Rrtmgp.H
+++ b/Source/Radiation/Rrtmgp.H
@@ -38,12 +38,13 @@ class Rrtmgp {
   public:
     // constructor
     Rrtmgp() = default;
-    explicit Rrtmgp(int ngases, char* gas_names[]);
     // deconstructor
     ~Rrtmgp() = default;
 
     // initialize and load gas property data for rrtmgp radiation
-    void initialize();
+    void initialize(int num_gas, std::vector<std::string> active_gas_names,
+                    const char* rrtmgp_coefficients_file_sw,
+                    const char* rrtmgp_coefficients_file_lw);
 
     // finalize/clean up
     void finalize();
@@ -117,7 +118,7 @@ class Rrtmgp {
 
    // number of gas for radiation model
    int ngas;
-   char** gas_names;
+   const char** gas_names;
 
    string1d active_gases;
 
diff --git a/Source/Radiation/data/cam5_ice_optics.nc b/Source/Radiation/data/cam5_ice_optics.nc
new file mode 100644
index 000000000..9d60ddc12
Binary files /dev/null and b/Source/Radiation/data/cam5_ice_optics.nc differ
diff --git a/Source/Radiation/data/cam5_liq_optics.nc b/Source/Radiation/data/cam5_liq_optics.nc
new file mode 100644
index 000000000..6dd6b077e
Binary files /dev/null and b/Source/Radiation/data/cam5_liq_optics.nc differ