Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Enable MODIS and MISR simulators from COSP #2838

Merged
merged 9 commits into from
Jun 3, 2024
56 changes: 34 additions & 22 deletions components/eamxx/src/physics/cosp/cosp_c2f.F90
Original file line number Diff line number Diff line change
Expand Up @@ -40,11 +40,11 @@ module cosp_c2f

! Local variables; control what runs and what does not
logical :: &
lsingle = .false., & ! True if using MMF_v3_single_moment CLOUDSAT microphysical scheme (default)
ldouble = .true., & ! True if using MMF_v3.5_two_moment CLOUDSAT microphysical scheme
lsingle = .false., & ! True if using MMF_v3_single_moment CLOUDSAT microphysical scheme (default)
ldouble = .true. , & ! True if using MMF_v3.5_two_moment CLOUDSAT microphysical scheme
lisccp = .true. , & ! Local on/off switch for simulators (used by initialization)
lmodis = .false., & !
lmisr = .false., & !
lmodis = .true. , & !
lmisr = .true. , & !
lcalipso = .false., & !
lgrLidar532 = .false., & !
latlid = .false., & !
Expand All @@ -64,7 +64,7 @@ module cosp_c2f
Lmeantbisccp = .false., & ! ISCCP mean all-sky 10.5micron brightness temperature
Lmeantbclrisccp = .false., & ! ISCCP mean clear-sky 10.5micron brightness temperature
Lalbisccp = .false., & ! ISCCP mean cloud albedo
LclMISR = .false., & ! MISR cloud fraction
LclMISR = .true. , & ! MISR cloud fraction
Lcltmodis = .false., & ! MODIS total cloud fraction
Lclwmodis = .false., & ! MODIS liquid cloud fraction
Lclimodis = .false., & ! MODIS ice cloud fraction
Expand All @@ -82,7 +82,7 @@ module cosp_c2f
Lpctmodis = .false., & ! MODIS cloud top pressure
Llwpmodis = .false., & ! MODIS cloud liquid water path
Liwpmodis = .false., & ! MODIS cloud ice water path
Lclmodis = .false., & ! MODIS cloud area fraction
Lclmodis = .true. , & ! MODIS cloud area fraction
Latb532 = .false., & ! CALIPSO attenuated total backscatter (532nm)
Latb532gr = .false., & ! GROUND LIDAR @ 532NM attenuated total backscatter (532nm)
Latb355 = .false., & ! ATLID attenuated total backscatter (355nm)
Expand Down Expand Up @@ -239,11 +239,9 @@ module cosp_c2f

subroutine cosp_c2f_init(npoints, ncolumns, nlevels) bind(c, name='cosp_c2f_init')
integer(kind=c_int), value, intent(in) :: npoints, ncolumns, nlevels
! Initialize/allocate COSP input and output derived types

! Number of vertical levels for Cloudsat/CALIPSO
nlvgrid = 40
call construct_cospIN(npoints,ncolumns,nlevels,cospIN)
call construct_cospstatein(npoints,nlevels,rttov_nchannels,cospstateIN)
call construct_cosp_outputs(npoints, ncolumns, nlevels, nlvgrid, rttov_nchannels, cospOUT)

! Initialize quickbeam_optics, also if two-moment radar microphysics scheme is wanted...
if (cloudsat_micro_scheme == 'MMF_v3.5_two_moment') then
Expand All @@ -264,19 +262,26 @@ subroutine cosp_c2f_init(npoints, ncolumns, nlevels) bind(c, name='cosp_c2f_init
cloudsat_radar_freq, cloudsat_k2, cloudsat_use_gas_abs, &
cloudsat_do_ray, isccp_topheight, isccp_topheight_direction, surface_radar, &
rcfg_cloudsat, use_vgrid, csat_vgrid, Nlvgrid, Nlevels, cloudsat_micro_scheme)

! Initialize/allocate COSP input and output derived types
call construct_cospIN(npoints,ncolumns,nlevels,cospIN)
call construct_cospstatein(npoints,nlevels,rttov_nchannels,cospstateIN)
call construct_cosp_outputs(npoints, ncolumns, nlevels, nlvgrid, rttov_nchannels, cospOUT)

end subroutine cosp_c2f_init

subroutine cosp_c2f_run(npoints, ncolumns, nlevels, ntau, nctp, &
emsfc_lw, sunlit, skt, T_mid, p_mid, p_int, qv, &
cldfrac, reff_qc, reff_qi, dtau067, dtau105, isccp_cldtot, isccp_ctptau &
subroutine cosp_c2f_run(npoints, ncolumns, nlevels, ntau, nctp, ncth, &
emsfc_lw, sunlit, skt, T_mid, p_mid, p_int, z_mid, qv, qc, qi, &
cldfrac, reff_qc, reff_qi, dtau067, dtau105, isccp_cldtot, isccp_ctptau, modis_ctptau, misr_cthtau &
) bind(C, name='cosp_c2f_run')
integer(kind=c_int), value, intent(in) :: npoints, ncolumns, nlevels, ntau, nctp
integer(kind=c_int), value, intent(in) :: npoints, ncolumns, nlevels, ntau, nctp, ncth
real(kind=c_double), value, intent(in) :: emsfc_lw
real(kind=c_double), intent(in), dimension(npoints) :: sunlit, skt
real(kind=c_double), intent(in), dimension(npoints,nlevels) :: T_mid, p_mid, qv, cldfrac, reff_qc, reff_qi, dtau067, dtau105
real(kind=c_double), intent(in), dimension(npoints,nlevels) :: T_mid, p_mid, z_mid, qv, qc, qi, cldfrac, reff_qc, reff_qi, dtau067, dtau105
real(kind=c_double), intent(in), dimension(npoints,nlevels+1) :: p_int
real(kind=c_double), intent(inout), dimension(npoints) :: isccp_cldtot
real(kind=c_double), intent(inout), dimension(npoints,ntau,nctp) :: isccp_ctptau
real(kind=c_double), intent(inout), dimension(npoints,ntau,nctp) :: isccp_ctptau, modis_ctptau
real(kind=c_double), intent(inout), dimension(npoints,ntau,ncth) :: misr_cthtau
! Takes normal arrays as input and populates COSP derived types
character(len=256),dimension(100) :: cosp_status
integer :: nptsperit
Expand Down Expand Up @@ -306,23 +311,25 @@ subroutine cosp_c2f_run(npoints, ncolumns, nlevels, ntau, nctp, &
cca(:npoints,:nlevels) = 0 ! Cloud fraction of convective clouds; not present or used in our model
dtau_c(:npoints,:nlevels) = 0 ! Optical depth of convective clouds; not present or used in our model
dem_c (:npoints,:nlevels) = 0 ! Emissivity of convective clouds; not present or used in our model
mr_lsliq(:npoints,:nlevels) = 0 ! Mixing ratio of cloud liquid; will be needed for radar/lidar
mr_ccliq(:npoints,:nlevels) = 0 ! Mixing ratio of cloud liquid for convective clouds; not present or used in our model
mr_lsice(:npoints,:nlevels) = 0 ! Mixing ratio of cloud ice; will be needed for radar/lidar
mr_lsliq(:npoints,:nlevels) = qc(:npoints,:nlevels) ! Mixing ratio of cloud liquid; will be needed for radar/lidar
mr_ccliq(:npoints,:nlevels) = 0 ! Mixing ratio of cloud liquid for convective clouds; not present or used in our model
mr_lsice(:npoints,:nlevels) = qi(:npoints,:nlevels) ! Mixing ratio of cloud ice; will be needed for radar/lidar
mr_ccice(:npoints,:nlevels) = 0 ! Mixing ratio of cloud ice for convective clouds; not present or used in our model
mr_lsrain(:npoints,:nlevels) = 0 ! Mixing ratio of rain; will be needed for radar/lidar
mr_ccrain(:npoints,:nlevels) = 0 ! Mixing ratio of rain for convective clouds; not present or used in our model
mr_lssnow(:npoints,:nlevels) = 0 ! Mixing ratio of snow; will be needed for radar/lidar
mr_ccsnow(:npoints,:nlevels) = 0 ! Mixing ratio of snow for convective clouds; will be needed for radar/lidar
mr_lsgrpl(:npoints,:nlevels) = 0 ! Mixing ratio of graupel; will be needed for radar/lidar
reff = 0 ! Effective radii; will be needed for MODIS
reff(:npoints,:nlevels,I_LSCLIQ) = 1e-6 * reff_qc(:npoints,:nlevels)
reff(:npoints,:nlevels,I_LSCICE) = 1e-6 * reff_qi(:npoints,:nlevels)

start_idx = 1
end_idx = npoints
! Translate arrays to derived types
cospIN%emsfc_lw = emsfc_lw
cospIN%rcfg_cloudsat = rcfg_cloudsat
! cospstateIN%hgt_matrix = zlev(start_idx:end_idx,Nlevels:1:-1) ! km
cospstateIN%hgt_matrix = z_mid(start_idx:end_idx,1:Nlevels) ! m
cospstateIN%sunlit = sunlit(start_idx:end_idx) ! 0-1
cospstateIN%skt = skt(start_idx:end_idx) ! K
! cospstateIN%surfelev = surfelev(start_idx:end_idx) ! m
Expand Down Expand Up @@ -351,6 +358,12 @@ subroutine cosp_c2f_run(npoints, ncolumns, nlevels, ntau, nctp, &
isccp_cldtot(:npoints) = cospOUT%isccp_totalcldarea(:npoints)
isccp_ctptau(:npoints,:,:) = cospOUT%isccp_fq(:npoints,:,:)

! Modis
modis_ctptau(:npoints,:,:) = cospOUT%modis_Optical_Thickness_vs_Cloud_Top_Pressure(:npoints,:,:)

! MISR
misr_cthtau(:npoints,:,:) = cospOUT%misr_fq(:npoints,:,:)

end subroutine cosp_c2f_run

subroutine cosp_c2f_final() bind(C, name='cosp_c2f_final')
Expand Down Expand Up @@ -455,8 +468,6 @@ end subroutine construct_cospstateIN
!
! This subroutine allocates output fields based on input logical flag switches.
! ######################################################################################
! TODO: This is WAY too many dummy arguments! These can just be defined at module scope I think
! and then initialized once at init
subroutine construct_cosp_outputs(Npoints,Ncolumns,Nlevels,Nlvgrid,Nchan,x)
! Inputs
integer,intent(in) :: &
Expand Down Expand Up @@ -1314,6 +1325,7 @@ subroutine subsample_and_optics(nPoints, nLevels, nColumns, nHydro, overlap, use

!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
! MODIS optics
! TODO: we can probably bypass this block if we take asym and ss_alb from radiation
!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if (Lmodis) then
allocate(MODIS_cloudWater(nPoints,nColumns,nLevels), &
Expand Down
31 changes: 21 additions & 10 deletions components/eamxx/src/physics/cosp/cosp_functions.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -4,11 +4,11 @@
using scream::Real;
extern "C" void cosp_c2f_init(int ncol, int nsubcol, int nlay);
extern "C" void cosp_c2f_final();
extern "C" void cosp_c2f_run(const int ncol, const int nsubcol, const int nlay, const int ntau, const int nctp,
extern "C" void cosp_c2f_run(const int ncol, const int nsubcol, const int nlay, const int ntau, const int nctp, const int ncth,
const Real emsfc_lw, const Real* sunlit, const Real* skt,
const Real* T_mid, const Real* p_mid, const Real* p_int, const Real* qv,
const Real* T_mid, const Real* p_mid, const Real* p_int, const Real* z_mid, const Real* qv, const Real* qc, const Real* qi,
const Real* cldfrac, const Real* reff_qc, const Real* reff_qi, const Real* dtau067, const Real* dtau105,
Real* isccp_cldtot, Real* isccp_ctptau);
Real* isccp_cldtot, Real* isccp_ctptau, Real* modis_ctptau, Real* misr_cthtau);

namespace scream {

Expand All @@ -30,28 +30,35 @@ namespace scream {
cosp_c2f_final();
};
inline void main(
const Int ncol, const Int nsubcol, const Int nlay, const Int ntau, const Int nctp, const Real emsfc_lw,
const Int ncol, const Int nsubcol, const Int nlay, const Int ntau, const Int nctp, const Int ncth, const Real emsfc_lw,
view_1d<const Real>& sunlit , view_1d<const Real>& skt,
view_2d<const Real>& T_mid , view_2d<const Real>& p_mid , view_2d<const Real>& p_int,
view_2d<const Real>& qv , view_2d<const Real>& cldfrac,
view_2d<const Real>& z_mid , view_2d<const Real>& qv , view_2d<const Real>& qc , view_2d<const Real>& qi,
view_2d<const Real>& cldfrac,
view_2d<const Real>& reff_qc, view_2d<const Real>& reff_qi,
view_2d<const Real>& dtau067, view_2d<const Real>& dtau105,
view_1d<Real>& isccp_cldtot , view_3d<Real>& isccp_ctptau) {
view_1d<Real>& isccp_cldtot , view_3d<Real>& isccp_ctptau, view_3d<Real>& modis_ctptau, view_3d<Real>& misr_cthtau) {

// Make host copies and permute data as needed
lview_host_2d
T_mid_h("T_mid_h", ncol, nlay), p_mid_h("p_mid_h", ncol, nlay), p_int_h("p_int_h", ncol, nlay+1),
qv_h("qv_h", ncol, nlay), cldfrac_h("cldfrac_h", ncol, nlay),
z_mid_h("z_mid_h", ncol, nlay), qv_h("qv_h", ncol, nlay), qc_h("qc_h", ncol, nlay), qi_h("qi_h", ncol, nlay),
cldfrac_h("cldfrac_h", ncol, nlay),
reff_qc_h("reff_qc_h", ncol, nlay), reff_qi_h("reff_qi_h", ncol, nlay),
dtau067_h("dtau_067_h", ncol, nlay), dtau105_h("dtau105_h", ncol, nlay);
lview_host_3d isccp_ctptau_h("isccp_ctptau_h", ncol, ntau, nctp);
lview_host_3d modis_ctptau_h("modis_ctptau_h", ncol, ntau, nctp);
lview_host_3d misr_cthtau_h("misr_cthtau_h", ncol, ntau, ncth);

// Copy to layoutLeft host views
for (int i = 0; i < ncol; i++) {
for (int j = 0; j < nlay; j++) {
T_mid_h(i,j) = T_mid(i,j);
p_mid_h(i,j) = p_mid(i,j);
z_mid_h(i,j) = z_mid(i,j);
qv_h(i,j) = qv(i,j);
qc_h(i,j) = qc(i,j);
qi_h(i,j) = qi(i,j);
cldfrac_h(i,j) = cldfrac(i,j);
reff_qc_h(i,j) = reff_qc(i,j);
reff_qi_h(i,j) = reff_qi(i,j);
Expand All @@ -68,17 +75,21 @@ namespace scream {
// Subsample here?

// Call COSP wrapper
cosp_c2f_run(ncol, nsubcol, nlay, ntau, nctp,
cosp_c2f_run(ncol, nsubcol, nlay, ntau, nctp, ncth,
emsfc_lw, sunlit.data(), skt.data(), T_mid_h.data(), p_mid_h.data(), p_int_h.data(),
qv_h.data(),
z_mid_h.data(), qv_h.data(), qc_h.data(), qi_h.data(),
cldfrac_h.data(), reff_qc_h.data(), reff_qi_h.data(), dtau067_h.data(), dtau105_h.data(),
isccp_cldtot.data(), isccp_ctptau_h.data());
isccp_cldtot.data(), isccp_ctptau_h.data(), modis_ctptau_h.data(), misr_cthtau_h.data());

// Copy outputs back to layoutRight views
for (int i = 0; i < ncol; i++) {
for (int j = 0; j < ntau; j++) {
for (int k = 0; k < nctp; k++) {
isccp_ctptau(i,j,k) = isccp_ctptau_h(i,j,k);
modis_ctptau(i,j,k) = modis_ctptau_h(i,j,k);
}
for (int k = 0; k < ncth; k++) {
misr_cthtau(i,j,k) = misr_cthtau_h(i,j,k);
}
}
}
Expand Down
Loading