LCOV - code coverage report
Current view: | top level - chemistry/utils - solar_irrad_data.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 91 | 119 | 76.5 %
Date: | 2024-12-17 17:57:11 |  | Functions: | 2 | 2 | 100.0 %

              Line data    Source code
           1             : !-----------------------------------------------------------------------
           2             : ! Solar irradiance / photon flux data
           3             : !-----------------------------------------------------------------------
           4             : module solar_irrad_data
           5             :   use shr_kind_mod,     only: r8 => shr_kind_r8
           6             :   use spmd_utils,       only: masterproc
           7             :   use cam_abortutils,   only: endrun
           8             :   use pio,              only: file_desc_t,pio_inq_dimid,pio_inq_varid,pio_inq_dimlen, pio_noerr,pio_internal_error,pio_bcast_error
           9             :   use pio,              only: pio_get_var, pio_seterrorhandling
          10             :   use cam_pio_utils,    only: cam_pio_openfile
          11             :   use cam_logfile,      only: iulog
          12             :   use infnan,           only: nan, assignment(=)
          13             :   use input_data_utils, only: time_coordinate
          14             :
          15             :   implicit none
          16             :
          17             :   save
          18             :   private
          19             :   public :: solar_irrad_init
          20             :   public :: solar_irrad_advance
          21             :
          22             :   integer,  public, protected :: nbins  ! number of wavelength samples of spectrum, wavelength endpoints
          23             :   real(r8), public, protected, allocatable :: we(:)
          24             :   real(r8), public, protected, allocatable :: sol_etf(:)
          25             :   real(r8), public, protected, allocatable :: ssi_ref(:)  ! a reference spectrum constructed from 3 solar cycles of data
          26             :   real(r8), public, protected, allocatable :: sol_irrad(:)
          27             :   real(r8), public, protected              :: sol_tsi = -1.0_r8
          28             :   real(r8), public, protected              :: ref_tsi
          29             :   logical,  public, protected :: do_spctrl_scaling = .false.
          30             :   logical,  public, protected :: has_spectrum = .false.
          31             :   logical,  public, protected :: has_ref_spectrum = .false.
          32             :
          33             :   type(file_desc_t) :: file_id
          34             :   integer :: ssi_vid
          35             :   integer :: tsi_vid
          36             :   integer :: ref_vid
          37             :   integer :: tsi_ref_vid
          38             :
          39             :   logical  :: initialized = .false.
          40             :   logical  :: has_tsi = .false.
          41             :   real(r8) :: itsi(2)
          42             :   real(r8), allocatable :: irradi(:,:)
          43             :   real(r8), allocatable :: irrad_fac(:)
          44             :   real(r8), allocatable :: etf_fac(:)
          45             :   real(r8), allocatable :: dellam(:)
          46             :
          47             :   logical  :: fixed_scon
          48             :
          49             :   type(time_coordinate) :: time_coord
          50             :
          51             : contains
          52             :
          53             : !-----------------------------------------------------------------------
          54             : !-----------------------------------------------------------------------
          55        1536 :   subroutine solar_irrad_init(filepath, fixed, fixed_ymd, fixed_tod, const_tsi, heatng_spctrl_scl )
          56             :
          57             :     use ioFileMod, only : getfil
          58             :     use physconst, only : c0, planck
          59             :
          60             :     !---------------------------------------------------------------
          61             :     ! arguments
          62             :     !---------------------------------------------------------------
          63             :     character(len=*), intent(in) :: filepath
          64             :     logical,  intent(in) :: fixed
          65             :     integer,  intent(in) :: fixed_ymd
          66             :     integer,  intent(in) :: fixed_tod
          67             :     real(r8), intent(in) :: const_tsi
          68             :     logical,  intent(in) :: heatng_spctrl_scl
          69             :
          70             :     !---------------------------------------------------------------
          71             :     ! local vars
          72             :     !---------------------------------------------------------------
          73             :     integer :: astat, dimid, vid
          74             :     character(len=256) :: filen
          75        1536 :     real(r8), allocatable :: lambda(:)
          76             :     integer :: i, wvl_vid
          77             :     real(r8), parameter :: c = c0     ! speed of light (m/s)
          78             :     real(r8), parameter :: h = planck ! Planck's constant (Joule sec)
          79             :     real(r8), parameter :: fac = 1._r8/(h*c)
          80             :     integer :: ierr
          81             :
          82        1536 :     has_spectrum = .false.
          83             :
          84        1536 :     if ( filepath.ne.'NONE' ) then
          85        1536 :        fixed_scon = .false.
          86             :     else
          87           0 :        fixed_scon = .true.
          88             :     endif
          89             :
          90        1536 :     if ( const_tsi>0._r8 ) then
          91           0 :        sol_tsi = const_tsi
          92             :     endif
          93        1536 :     ref_tsi = nan
          94             :
          95        1536 :     if ( fixed_scon ) return
          96             :
          97             :     call time_coord%initialize( filepath, fixed=fixed, fixed_ymd=fixed_ymd, fixed_tod=fixed_tod, &
          98        1536 :                                 force_time_interp=.true., try_dates=.true. )
          99             :
         100        1536 :     call getfil( filepath, filen, 0 )
         101        1536 :     call cam_pio_openfile( file_id, filen, 0 )
         102        1536 :     if(masterproc)   write(iulog,*)'solar_data_init: data file = ',trim(filen)
         103        1536 :     call pio_seterrorhandling(file_id, pio_bcast_error)
         104        1536 :     ierr = pio_inq_varid( file_id, 'ssi', ssi_vid )
         105        1536 :     has_spectrum = ierr==PIO_NOERR
         106             :
         107        1536 :     ierr = pio_inq_varid( file_id, 'tsi', tsi_vid )
         108        1536 :     has_tsi = ierr==PIO_NOERR .and. const_tsi<0._r8
         109             :
         110        1536 :     ierr = pio_inq_varid( file_id, 'ssi_ref', ref_vid )
         111        1536 :     has_ref_spectrum = ierr==PIO_NOERR
         112        1536 :     call pio_seterrorhandling(file_id, pio_internal_error)
         113             :
         114        1536 :     if ( has_spectrum ) then
         115        1536 :        call pio_seterrorhandling(file_id, pio_bcast_error)
         116        1536 :        ierr = pio_inq_varid( file_id, 'wavelength', wvl_vid )
         117        1536 :        call pio_seterrorhandling(file_id, pio_internal_error)
         118             :
         119        1536 :        if ( ierr==PIO_NOERR ) then
         120        1536 :           ierr = pio_inq_dimid( file_id, 'wavelength', dimid )
         121             :        else ! for backwards compatibility
         122           0 :           ierr = pio_inq_varid( file_id, 'wvl', wvl_vid  )
         123           0 :           ierr = pio_inq_dimid( file_id, 'wvl', dimid )
         124             :        endif
         125        1536 :        ierr = pio_inq_dimlen( file_id, dimid, nbins )
         126        1536 :        if ( has_ref_spectrum ) then
         127           0 :           ierr = pio_inq_varid( file_id, 'tsi_ref', tsi_ref_vid )
         128             :        endif
         129             :     endif
         130             :
         131        1536 :     do_spctrl_scaling = has_spectrum .and. heatng_spctrl_scl
         132             :
         133        4608 :     allocate(lambda(nbins), stat=astat )
         134        1536 :     if( astat /= 0 ) then
         135           0 :        write(iulog,*) 'solar_data_init: failed to allocate lambda; error = ',astat
         136           0 :        call endrun('solar_data_init')
         137             :     end if
         138        4608 :     allocate(dellam(nbins), stat=astat )
         139        1536 :     if( astat /= 0 ) then
         140           0 :        write(iulog,*) 'solar_data_init: failed to allocate dellam; error = ',astat
         141           0 :        call endrun('solar_data_init')
         142             :     end if
         143        4608 :     allocate(irrad_fac(nbins), stat=astat )
         144        1536 :     if( astat /= 0 ) then
         145           0 :        write(iulog,*) 'solar_data_init: failed to allocate irrad_fac; error = ',astat
         146           0 :        call endrun('solar_data_init')
         147             :     end if
         148        4608 :     allocate(etf_fac(nbins), stat=astat )
         149        1536 :     if( astat /= 0 ) then
         150           0 :        write(iulog,*) 'solar_data_init: failed to allocate etf_fac; error = ',astat
         151           0 :        call endrun('solar_data_init')
         152             :     end if
         153        4608 :     allocate(sol_irrad(nbins), stat=astat )
         154        1536 :     if( astat /= 0 ) then
         155           0 :        write(iulog,*) 'solar_data_init: failed to allocate sol_irrad; error = ',astat
         156           0 :        call endrun('solar_data_init')
         157             :     end if
         158        4608 :     allocate(ssi_ref(nbins), stat=astat )
         159        1536 :     if( astat /= 0 ) then
         160           0 :        write(iulog,*) 'solar_data_init: failed to allocate ssi_ref; error = ',astat
         161           0 :        call endrun('solar_data_init')
         162             :     end if
         163        1536 :     ssi_ref(:) = nan
         164             :
         165        1536 :     if (has_spectrum) then
         166        1536 :        ierr = pio_get_var( file_id, wvl_vid, lambda )
         167        1536 :        ierr = pio_inq_varid( file_id, 'band_width', vid  )
         168        1536 :        ierr = pio_get_var( file_id, vid, dellam )
         169             :     endif
         170             :
         171        1536 :     if(masterproc)   write(iulog,*)'solar_data_init: has_ref_spectrum',has_ref_spectrum
         172        1536 :     if ( has_ref_spectrum ) then
         173           0 :        ierr = pio_inq_varid( file_id, 'ssi_ref', vid  )
         174           0 :        ierr = pio_get_var( file_id, vid, ssi_ref )
         175           0 :        ierr = pio_get_var( file_id, tsi_ref_vid, ref_tsi )
         176             :     endif
         177             :
         178        1536 :     if ( has_spectrum ) then
         179        4608 :        allocate(sol_etf(nbins), stat=astat )
         180        1536 :        if( astat /= 0 ) then
         181           0 :           write(iulog,*) 'solar_data_init: failed to allocate sol_etf; error = ',astat
         182           0 :           call endrun('solar_data_init')
         183             :        end if
         184        4608 :        allocate(irradi(nbins,2), stat=astat )
         185        1536 :        if( astat /= 0 ) then
         186           0 :           write(iulog,*) 'solar_data_init: failed to allocate irradi; error = ',astat
         187           0 :           call endrun('solar_data_init')
         188             :        end if
         189             :
         190        4608 :        allocate(we(nbins+1), stat=astat )
         191        1536 :        if( astat /= 0 ) then
         192           0 :           write(iulog,*) 'solar_data_init: failed to allocate we; error = ',astat
         193           0 :           call endrun('solar_data_init')
         194             :        end if
         195             :
         196     5976576 :        we(:nbins)  = lambda(:nbins) - 0.5_r8*dellam(:nbins)
         197        1536 :        we(nbins+1) = lambda(nbins)  + 0.5_r8*dellam(nbins)
         198     5976576 :        do i = 1,nbins
         199     5975040 :           irrad_fac(i) = 1.e-3_r8                ! mW/m2/nm --> W/m2/nm
         200     5976576 :           etf_fac(i)   = 1.e-16_r8*lambda(i)*fac ! mW/m2/nm --> photons/cm2/sec/nm
         201             :        enddo
         202        1536 :        if(has_ref_spectrum) then
         203           0 :           ssi_ref = ssi_ref * 1.e-3_r8        ! mW/m2/nm --> W/m2/nm
         204             :        endif
         205             :     endif
         206             :
         207        1536 :     deallocate(lambda)
         208        1536 :     deallocate(dellam)
         209             :
         210             :     ! need to force data loading when the model starts at a time =/ 00:00:00.000
         211             :     ! -- may occur in restarts also
         212        1536 :     call solar_irrad_advance()
         213        1536 :     initialized = .true.
         214             :
         215        1536 :   end subroutine solar_irrad_init
         216             :
         217             : !-----------------------------------------------------------------------
         218             : ! Reads in the ETF data for the current date.
         219             : !-----------------------------------------------------------------------
         220      372480 :   subroutine solar_irrad_advance( )
         221             :
         222             :     integer  :: i, index, nt
         223             :     integer  :: offset(2), count(2)
         224             :     logical  :: read_data
         225      372480 :     real(r8) :: data(nbins)
         226             :     integer  :: ierr
         227             :     real(r8) :: delt
         228             :
         229      372480 :     if ( fixed_scon ) return
         230      372480 :     if ( time_coord%fixed .and. initialized ) return
         231             :
         232      372480 :     index = -1
         233             :
         234      372480 :     read_data = time_coord%read_more() .or. .not.initialized
         235      372480 :     call time_coord%advance()
         236             :
         237      372480 :     if ( read_data ) then
         238        9216 :        nt = 2
         239        9216 :        index = time_coord%indxs(1)
         240             :
         241             :        ! get the surrounding time slices
         242       27648 :        offset = (/ 1, index /)
         243       27648 :        count =  (/ nbins, nt /)
         244             :
         245        9216 :        if (has_spectrum) then
         246        9216 :           ierr = pio_get_var( file_id, ssi_vid, offset, count, irradi )
         247             :        endif
         248        9216 :        if (has_tsi .and. (.not.do_spctrl_scaling)) then
         249       27648 :           ierr = pio_get_var( file_id, tsi_vid, (/index/), (/nt/), itsi )
         250       27648 :           if ( any(itsi(:nt) < 0._r8) ) then
         251           0 :              call endrun( 'solar_data_advance: invalid or missing tsi data  ' )
         252             :           endif
         253             :        endif
         254             :     endif
         255             :
         256      372480 :     delt = time_coord%wghts(2)
         257             :
         258      372480 :     if (has_spectrum) then
         259  1449319680 :        data(:) = irradi(:,1) + delt*( irradi(:,2) - irradi(:,1) )
         260             :
         261  1449319680 :        do i = 1,nbins
         262  1448947200 :           sol_irrad(i) = data(i)*irrad_fac(i) ! W/m2/nm
         263  1449319680 :           sol_etf(i)   = data(i)*etf_fac(i)   ! photons/cm2/sec/nm
         264             :        enddo
         265             :     endif
         266      372480 :     if (has_tsi .and. (.not.do_spctrl_scaling)) then
         267      372480 :        sol_tsi = itsi(1) + delt*( itsi(2) - itsi(1) )
         268             :     endif
         269             :
         270             :   end subroutine solar_irrad_advance
         271             :
         272             : end module solar_irrad_data

Generated by: LCOV version 1.14