Skip to content

Commit

Permalink
Enable modis_ctptau calculation and output
Browse files Browse the repository at this point in the history
  • Loading branch information
@brhillman
brhillman committed May 18, 2024
1 parent 8ff32a8 commit 4c11ca5
Show file tree
Hide file tree
Showing 3 changed files with 28 additions and 16 deletions.
25 changes: 15 additions & 10 deletions components/eamxx/src/physics/cosp/cosp_c2f.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,7 @@ module cosp_c2f
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., & !
lmodis = .true., & !
lmisr = .false., & !
lcalipso = .false., & !
lgrLidar532 = .false., & !
Expand Down Expand Up @@ -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,25 @@ 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 &
cldfrac, reff_qc, reff_qi, dtau067, dtau105, isccp_cldtot, isccp_ctptau, modis_ctptau &
) bind(C, name='cosp_c2f_run')
integer(kind=c_int), value, intent(in) :: npoints, ncolumns, nlevels, ntau, nctp
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+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
! Takes normal arrays as input and populates COSP derived types
character(len=256),dimension(100) :: cosp_status
integer :: nptsperit
Expand Down Expand Up @@ -351,6 +355,9 @@ 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,:,:)

end subroutine cosp_c2f_run

subroutine cosp_c2f_final() bind(C, name='cosp_c2f_final')
Expand Down Expand Up @@ -455,8 +462,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
8 changes: 5 additions & 3 deletions components/eamxx/src/physics/cosp/cosp_functions.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@ extern "C" void cosp_c2f_run(const int ncol, const int nsubcol, const int nlay,
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* 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);

namespace scream {

Expand Down Expand Up @@ -36,7 +36,7 @@ namespace scream {
view_2d<const Real>& qv , 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) {

// Make host copies and permute data as needed
lview_host_2d
Expand All @@ -45,6 +45,7 @@ namespace scream {
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);

// Copy to layoutLeft host views
for (int i = 0; i < ncol; i++) {
Expand Down Expand Up @@ -72,13 +73,14 @@ namespace scream {
emsfc_lw, sunlit.data(), skt.data(), T_mid_h.data(), p_mid_h.data(), p_int_h.data(),
qv_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());

// 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);
}
}
}
Expand Down
11 changes: 8 additions & 3 deletions components/eamxx/src/physics/cosp/eamxx_cosp.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -85,6 +85,7 @@ void Cosp::set_grids(const std::shared_ptr<const GridsManager> grids_manager)
// Set of fields used strictly as output
add_field<Computed>("isccp_cldtot", scalar2d, percent, grid_name);
add_field<Computed>("isccp_ctptau", scalar4d_ctptau, percent, grid_name, 1);
add_field<Computed>("modis_ctptau", scalar4d_ctptau, percent, grid_name, 1);
add_field<Computed>("isccp_mask" , scalar2d, nondim, grid_name);
}

Expand All @@ -98,8 +99,8 @@ void Cosp::initialize_impl (const RunType /* run_type */)
// Add note to output files about processing ISCCP fields that are only valid during
// daytime. This can go away once I/O can handle masked time averages.
using stratts_t = std::map<std::string,std::string>;
std::list<std::string> vnames = {"isccp_cldtot", "isccp_ctptau"};
for (const auto field_name : {"isccp_cldtot", "isccp_ctptau"}) {
std::list<std::string> vnames = {"isccp_cldtot", "isccp_ctptau", "modis_ctptau"};
for (const auto field_name : {"isccp_cldtot", "isccp_ctptau", "modis_ctptau"}) {
auto& f = get_field_out(field_name);
auto& atts = f.get_header().get_extra_data<stratts_t>("io: string attributes");
atts["note"] = "Night values are zero; divide by isccp_mask to get daytime mean";
Expand Down Expand Up @@ -162,6 +163,7 @@ void Cosp::run_impl (const double dt)
auto dtau105 = get_field_in("dtau105").get_view<const Real**, Host>();
auto isccp_cldtot = get_field_out("isccp_cldtot").get_view<Real*, Host>();
auto isccp_ctptau = get_field_out("isccp_ctptau").get_view<Real***, Host>();
auto modis_ctptau = get_field_out("modis_ctptau").get_view<Real***, Host>();
auto isccp_mask = get_field_out("isccp_mask" ).get_view<Real*, Host>(); // Copy of sunlit flag with COSP frequency for proper averaging

// Call COSP wrapper routines
Expand All @@ -172,7 +174,7 @@ void Cosp::run_impl (const double dt)
m_num_cols, m_num_subcols, m_num_levs, m_num_isccptau, m_num_isccpctp,
emsfc_lw, sunlit, skt, T_mid, p_mid, p_int, qv,
cldfrac, reff_qc, reff_qi, dtau067, dtau105,
isccp_cldtot, isccp_ctptau
isccp_cldtot, isccp_ctptau, modis_ctptau
);
// Remask night values to ZERO since our I/O does not know how to handle masked/missing values
// in temporal averages; this is all host data, so we can just use host loops like its the 1980s
Expand All @@ -182,6 +184,7 @@ void Cosp::run_impl (const double dt)
for (int j = 0; j < m_num_isccptau; j++) {
for (int k = 0; k < m_num_isccpctp; k++) {
isccp_ctptau(i,j,k) = 0;
modis_ctptau(i,j,k) = 0;
}
}
}
Expand All @@ -198,10 +201,12 @@ void Cosp::run_impl (const double dt)
// TODO: mask this when/if the AD ever supports masked averages
Kokkos::deep_copy(isccp_cldtot, 0.0);
Kokkos::deep_copy(isccp_ctptau, 0.0);
Kokkos::deep_copy(modis_ctptau, 0.0);
Kokkos::deep_copy(isccp_mask , 0.0);
}
get_field_out("isccp_cldtot").sync_to_dev();
get_field_out("isccp_ctptau").sync_to_dev();
get_field_out("modis_ctptau").sync_to_dev();
get_field_out("isccp_mask" ).sync_to_dev();
}

Expand Down

0 comments on commit 4c11ca5

Please sign in to comment.