|
| 1 | +!------------------------------------------------------------------------------- |
| 2 | +! This module uses the solar irradiance data |
| 3 | +! to provide a spectral scaling factor |
| 4 | +! to approximate the spectral distribution of irradiance |
| 5 | +! when the radiation scheme might use a different solar source function |
| 6 | +!------------------------------------------------------------------------------- |
| 7 | +module rad_solar_var |
| 8 | + |
| 9 | + use shr_kind_mod , only : r8 => shr_kind_r8 |
| 10 | + use radconstants, only : nswbands, get_sw_spectral_boundaries, band2gpt_sw |
| 11 | + use solar_irrad_data, only : sol_irrad, we, nbins, has_spectrum, sol_tsi |
| 12 | + use solar_irrad_data, only : do_spctrl_scaling |
| 13 | + use cam_abortutils, only : endrun |
| 14 | + use error_messages, only : alloc_err |
| 15 | + |
| 16 | + implicit none |
| 17 | + save |
| 18 | + |
| 19 | + private |
| 20 | + public :: rad_solar_var_init |
| 21 | + public :: get_variability |
| 22 | + |
| 23 | + real(r8), allocatable :: irrad(:) ! solar irradiance at model timestep in each band |
| 24 | + |
| 25 | + real(r8), allocatable :: radbinmax(:) |
| 26 | + real(r8), allocatable :: radbinmin(:) |
| 27 | + |
| 28 | +!------------------------------------------------------------------------------- |
| 29 | +contains |
| 30 | +!------------------------------------------------------------------------------- |
| 31 | + |
| 32 | + subroutine rad_solar_var_init( ) |
| 33 | + |
| 34 | + integer :: ierr |
| 35 | + integer :: radmax_loc |
| 36 | + |
| 37 | + if ( do_spctrl_scaling ) then |
| 38 | + |
| 39 | + if ( .not.has_spectrum ) then |
| 40 | + call endrun('rad_solar_var_init: solar input file must have irradiance spectrum') |
| 41 | + endif |
| 42 | + |
| 43 | + allocate (radbinmax(nswbands),stat=ierr) |
| 44 | + if (ierr /= 0) then |
| 45 | + call endrun('rad_solar_var_init: Error allocating space for radbinmax') |
| 46 | + end if |
| 47 | + |
| 48 | + allocate (radbinmin(nswbands),stat=ierr) |
| 49 | + if (ierr /= 0) then |
| 50 | + call endrun('rad_solar_var_init: Error allocating space for radbinmin') |
| 51 | + end if |
| 52 | + |
| 53 | + allocate (irrad(nswbands), stat=ierr) |
| 54 | + if (ierr /= 0) then |
| 55 | + call endrun('rad_solar_var_init: Error allocating space for irrad') |
| 56 | + end if |
| 57 | + |
| 58 | + call get_sw_spectral_boundaries(radbinmin, radbinmax, 'nm') |
| 59 | + |
| 60 | + ! Make sure that the far-IR is included, even if radiation grid does not |
| 61 | + ! extend that far down. 10^5 nm corresponds to a wavenumber of |
| 62 | + ! 100 cm^-1. |
| 63 | + radmax_loc = maxloc(radbinmax,1) |
| 64 | + radbinmax(radmax_loc) = max(100000._r8,radbinmax(radmax_loc)) |
| 65 | + |
| 66 | + endif |
| 67 | + |
| 68 | + end subroutine rad_solar_var_init |
| 69 | + |
| 70 | +!------------------------------------------------------------------------------- |
| 71 | +!------------------------------------------------------------------------------- |
| 72 | + |
| 73 | + subroutine get_variability(toa_flux, sfac) |
| 74 | + |
| 75 | + ! Arguments |
| 76 | + real(r8), intent(in) :: toa_flux(:,:) ! TOA flux to be scaled (columns,gpts) |
| 77 | + real(r8), intent(out) :: sfac(:,:) ! scaling factors (columns,gpts) |
| 78 | + |
| 79 | + ! Local variables |
| 80 | + integer :: i, j, istat, gpt_start, gpt_end, ncols |
| 81 | + real(r8), allocatable :: scale(:) |
| 82 | + character(len=*), parameter :: sub = 'get_variability' |
| 83 | + |
| 84 | + if (do_spctrl_scaling) then |
| 85 | + |
| 86 | + ! Determine target irradiance for each band |
| 87 | + call integrate_spectrum(nbins, nswbands, we, radbinmin, radbinmax, sol_irrad, irrad) |
| 88 | + |
| 89 | + ncols = size(toa_flux, 1) |
| 90 | + allocate(scale(ncols), stat=istat) |
| 91 | + call alloc_err(istat, sub, 'scale', ncols) |
| 92 | + |
| 93 | + do i = 1, nswbands |
| 94 | + gpt_start = band2gpt_sw(1,i) |
| 95 | + gpt_end = band2gpt_sw(2,i) |
| 96 | + scale = spread(irrad(i), 1, ncols) / sum(toa_flux(:, gpt_start:gpt_end), dim=2) |
| 97 | + do j = gpt_start, gpt_end |
| 98 | + sfac(:,j) = scale |
| 99 | + end do |
| 100 | + end do |
| 101 | + |
| 102 | + else |
| 103 | + sfac(:,:) = sol_tsi / spread(sum(toa_flux, 2), 2, size(toa_flux, 2)) |
| 104 | + end if |
| 105 | + end subroutine get_variability |
| 106 | + |
| 107 | + |
| 108 | +!------------------------------------------------------------------------------- |
| 109 | +! private method......... |
| 110 | +!------------------------------------------------------------------------------- |
| 111 | + |
| 112 | + subroutine integrate_spectrum( nsrc, ntrg, src_x, min_trg, max_trg, src, trg ) |
| 113 | + |
| 114 | + use mo_util, only : rebin |
| 115 | + |
| 116 | + implicit none |
| 117 | + |
| 118 | + !--------------------------------------------------------------- |
| 119 | + ! ... dummy arguments |
| 120 | + !--------------------------------------------------------------- |
| 121 | + integer, intent(in) :: nsrc ! dimension source array |
| 122 | + integer, intent(in) :: ntrg ! dimension target array |
| 123 | + real(r8), intent(in) :: src_x(nsrc+1) ! source coordinates |
| 124 | + real(r8), intent(in) :: max_trg(ntrg) ! target coordinates |
| 125 | + real(r8), intent(in) :: min_trg(ntrg) ! target coordinates |
| 126 | + real(r8), intent(in) :: src(nsrc) ! source array |
| 127 | + real(r8), intent(out) :: trg(ntrg) ! target array |
| 128 | + |
| 129 | + !--------------------------------------------------------------- |
| 130 | + ! ... local variables |
| 131 | + !--------------------------------------------------------------- |
| 132 | + real(r8) :: trg_x(2), targ(1) ! target coordinates |
| 133 | + integer :: i |
| 134 | + |
| 135 | + do i = 1, ntrg |
| 136 | + |
| 137 | + trg_x(1) = min_trg(i) |
| 138 | + trg_x(2) = max_trg(i) |
| 139 | + |
| 140 | + call rebin( nsrc, 1, src_x, trg_x, src(1:nsrc), targ(:) ) |
| 141 | + ! W/m2/nm --> W/m2 |
| 142 | + trg( i ) = targ(1)*(trg_x(2)-trg_x(1)) |
| 143 | + |
| 144 | + enddo |
| 145 | + |
| 146 | + |
| 147 | + end subroutine integrate_spectrum |
| 148 | + |
| 149 | +end module rad_solar_var |
0 commit comments