LCOV - code coverage report
Current view: | top level - chemistry/mozart - mee_fluxes.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 0 | 101 | 0.0 %
Date: | 2025-01-13 21:54:50 |  | Functions: | 0 | 6 | 0.0 %

              Line data    Source code
           1             : !--------------------------------------------------------------------------------
           2             : ! Provids electron fluxes read from input data set
           3             : !--------------------------------------------------------------------------------
           4             : module mee_fluxes
           5             :   use shr_kind_mod,   only : r8 => shr_kind_r8, cl=> shr_kind_cl
           6             :   use spmd_utils,     only : masterproc
           7             :   use cam_logfile,    only : iulog
           8             :   use cam_abortutils, only : endrun
           9             :   use input_data_utils, only : time_coordinate
          10             :   use pio, only : file_desc_t, var_desc_t, pio_get_var, pio_inq_varid
          11             :   use pio, only : PIO_NOWRITE, pio_inq_dimid, pio_inq_dimlen
          12             :   use infnan, only: isnan
          13             :
          14             :   implicit none
          15             :
          16             :   private
          17             :   public :: mee_fluxes_readnl
          18             :   public :: mee_fluxes_init
          19             :   public :: mee_fluxes_final
          20             :   public :: mee_fluxes_adv     ! read and time interpolate fluxes
          21             :   public :: mee_fluxes_extract ! interpolate flux to column L-shell
          22             :   public :: mee_fluxes_active  ! true when input flux file is specified
          23             :   public :: mee_fluxes_denergy ! energy bin widths
          24             :   public :: mee_fluxes_energy  ! center of each energy bin
          25             :   public :: mee_fluxes_nenergy ! number of energy bins
          26             :
          27             :   real(r8),protected, pointer :: mee_fluxes_denergy(:) => null()
          28             :   real(r8),protected, pointer :: mee_fluxes_energy(:) => null()
          29             :   integer, protected :: mee_fluxes_nenergy
          30             :   logical, protected :: mee_fluxes_active = .false.
          31             :
          32             :   real(r8), allocatable :: lshell(:)
          33             :   real(r8), allocatable :: indata(:,:,:)
          34             :   real(r8), allocatable :: influx(:,:)
          35             :   logical , allocatable :: valflx(:,:)
          36             :
          37             :   character(len=cl) :: mee_fluxes_filepath = 'NONE'
          38             :   logical :: mee_fluxes_fillin = .false.
          39             :
          40             :   type(time_coordinate) :: time_coord
          41             :   integer :: nlshells
          42             :
          43             :   type(file_desc_t) :: file_id
          44             :   type(var_desc_t) :: flux_var_id
          45             :
          46             : contains
          47             :
          48             :   !-----------------------------------------------------------------------------
          49             :   ! read namelist options
          50             :   !-----------------------------------------------------------------------------
          51           0 :   subroutine mee_fluxes_readnl(nlfile)
          52             :
          53             :     use namelist_utils, only: find_group_name
          54             :     use spmd_utils,     only: mpicom, mpi_character, mpi_logical, masterprocid
          55             :
          56             :     character(len=*), intent(in) :: nlfile  ! filepath for file containing namelist input
          57             :
          58             :     ! Local variables
          59             :     integer :: unitn, ierr
          60             :     character(len=*), parameter :: subname = 'mee_fluxes_readnl'
          61             :
          62             :     namelist /mee_fluxes_opts/ mee_fluxes_filepath, mee_fluxes_fillin
          63             :
          64             :     ! Read namelist
          65           0 :     if (masterproc) then
          66           0 :        open( newunit=unitn, file=trim(nlfile), status='old' )
          67           0 :        call find_group_name(unitn, 'mee_fluxes_opts', status=ierr)
          68           0 :        if (ierr == 0) then
          69           0 :           read(unitn, mee_fluxes_opts, iostat=ierr)
          70           0 :           if (ierr /= 0) then
          71           0 :              call endrun(subname // ':: ERROR reading namelist')
          72             :           end if
          73             :        end if
          74           0 :        close(unitn)
          75             :     end if
          76             :
          77             :     ! Broadcast namelist variables
          78           0 :     call mpi_bcast(mee_fluxes_filepath, len(mee_fluxes_filepath), mpi_character, masterprocid, mpicom, ierr)
          79           0 :     call mpi_bcast(mee_fluxes_fillin, 1, mpi_logical, masterprocid, mpicom, ierr)
          80             :
          81           0 :     mee_fluxes_active = mee_fluxes_filepath /= 'NONE'
          82             :
          83           0 :     if ( masterproc ) then
          84           0 :        if ( mee_fluxes_active ) then
          85           0 :           write(iulog,*) subname//':: Input electron fluxes filepath: '//trim(mee_fluxes_filepath)
          86           0 :           write(iulog,*) subname//':: Fill in missing fluxes with vdk-derived fluxes: ', mee_fluxes_fillin
          87             :        else
          88           0 :           write(iulog,*) subname//':: Electron fluxes are not prescribed'
          89             :        end if
          90             :     end if
          91             :
          92           0 :   end subroutine mee_fluxes_readnl
          93             :
          94             :   !-----------------------------------------------------------------------------
          95             :   ! intialize -- allocate memory and read coordinate data
          96             :   !-----------------------------------------------------------------------------
          97           0 :   subroutine mee_fluxes_init()
          98             :     use cam_pio_utils,  only : cam_pio_openfile
          99             :     use ioFileMod,      only : getfil
         100             :
         101             :     character(len=cl) :: filen
         102             :     integer :: ierr, dimid, varid
         103           0 :     real(r8), allocatable :: logdelta(:)
         104             :     character(len=*), parameter :: subname = 'mee_fluxes_init: '
         105             :
         106           0 :     if (.not.mee_fluxes_active) return
         107             :
         108           0 :     call time_coord%initialize( mee_fluxes_filepath, force_time_interp=.true. )
         109             :
         110           0 :     call getfil( mee_fluxes_filepath, filen, 0 )
         111           0 :     call cam_pio_openfile( file_id, filen, PIO_NOWRITE )
         112             :
         113           0 :     ierr = pio_inq_dimid(file_id, 'energy', dimid)
         114           0 :     ierr = pio_inq_dimlen(file_id, dimid, mee_fluxes_nenergy )
         115             :
         116           0 :     ierr = pio_inq_dimid(file_id, 'lshell', dimid)
         117           0 :     ierr = pio_inq_dimlen(file_id, dimid, nlshells )
         118             :
         119           0 :     ierr = pio_inq_varid(file_id, 'RBSP_flux_scaled', flux_var_id)
         120             :
         121           0 :     allocate( indata( mee_fluxes_nenergy, nlshells, 2 ), stat=ierr)
         122           0 :     if (ierr/=0) call endrun(subname//'not able to allocate indata')
         123           0 :     allocate( influx( mee_fluxes_nenergy, nlshells ), stat=ierr )
         124           0 :     if (ierr/=0) call endrun(subname//'not able to allocate influx')
         125           0 :     allocate( valflx( mee_fluxes_nenergy, nlshells ), stat=ierr )
         126           0 :     if (ierr/=0) call endrun(subname//'not able to allocate valflx')
         127           0 :     allocate( mee_fluxes_energy( mee_fluxes_nenergy ), stat=ierr )
         128           0 :     if (ierr/=0) call endrun(subname//'not able to allocate mee_fluxes_energy')
         129           0 :     allocate( mee_fluxes_denergy( mee_fluxes_nenergy ), stat=ierr )
         130           0 :     if (ierr/=0) call endrun(subname//'not able to allocate mee_fluxes_denergy')
         131           0 :     allocate( logdelta( mee_fluxes_nenergy ), stat=ierr )
         132           0 :     if (ierr/=0) call endrun(subname//'not able to allocate logdelta')
         133           0 :     allocate( lshell( nlshells ), stat=ierr )
         134           0 :     if (ierr/=0) call endrun(subname//'not able to allocate lshell')
         135             :
         136             :
         137           0 :     ierr = pio_inq_varid(file_id, 'energy', varid)
         138           0 :     ierr = pio_get_var( file_id, varid, mee_fluxes_energy)
         139             :
         140           0 :     ierr = pio_inq_varid(file_id, 'lshell', varid)
         141           0 :     ierr = pio_get_var( file_id, varid, lshell)
         142             :
         143           0 :     logdelta(2:) =  log(mee_fluxes_energy(2:mee_fluxes_nenergy))-log(mee_fluxes_energy(1:mee_fluxes_nenergy-1))
         144           0 :     logdelta(1) = logdelta(2)
         145           0 :     mee_fluxes_denergy(:) = exp( log(mee_fluxes_energy(:)) + 0.5_r8*logdelta(:) ) &
         146           0 :                           - exp( log(mee_fluxes_energy(:)) - 0.5_r8*logdelta(:) )
         147             :
         148           0 :     deallocate(logdelta)
         149             :
         150           0 :     call read_fluxes()
         151             :
         152           0 :   end subroutine mee_fluxes_init
         153             :
         154             :   !-----------------------------------------------------------------------------
         155             :   !-----------------------------------------------------------------------------
         156           0 :   subroutine mee_fluxes_final()
         157           0 :     use cam_pio_utils, only : cam_pio_closefile
         158             :
         159           0 :     if (.not.mee_fluxes_active) return
         160             :
         161           0 :     call cam_pio_closefile(file_id)
         162             :
         163           0 :     deallocate(indata)
         164           0 :     deallocate(influx)
         165           0 :     deallocate(valflx)
         166           0 :     deallocate(lshell)
         167             :
         168           0 :     deallocate(mee_fluxes_energy)
         169           0 :     deallocate(mee_fluxes_denergy)
         170             :
         171           0 :     nullify(mee_fluxes_energy)
         172             :     nullify(mee_fluxes_denergy)
         173             :
         174           0 :   end subroutine mee_fluxes_final
         175             :
         176             :   !-----------------------------------------------------------------------------
         177             :   ! time interpolate the input fluxes
         178             :   ! reads data as needed
         179             :   !-----------------------------------------------------------------------------
         180           0 :   subroutine mee_fluxes_adv
         181             :
         182           0 :     if (.not.mee_fluxes_active) return
         183             :
         184           0 :     call time_coord%advance()
         185             :
         186           0 :     if ( time_coord%read_more() ) then
         187           0 :        call read_fluxes( )
         188             :     endif
         189             :
         190           0 :     influx(:,:) = 0._r8
         191             :
         192           0 :     valflx(:,:) = (.not.isnan(indata(:,:,1))) .and. (.not.isnan(indata(:,:,2)))
         193           0 :     where ( valflx(:,:) )
         194           0 :        influx(:,:) = time_coord%wghts(1)*indata(:,:,1) + time_coord%wghts(2)*indata(:,:,2)
         195             :     end where
         196             :
         197           0 :     if (any(isnan(influx))) then
         198           0 :        call endrun('mee_fluxes_adv -- influx has NaNs')
         199             :     end if
         200             :
         201           0 :   end subroutine mee_fluxes_adv
         202             :
         203             :   !-----------------------------------------------------------------------------
         204             :   ! linear interpolate fluxes in L-shell where the fluxes are valid
         205             :   !-----------------------------------------------------------------------------
         206           0 :   subroutine mee_fluxes_extract( l_shell, fluxes, valid )
         207             :
         208             :     real(r8), intent(in) :: l_shell
         209             :     real(r8), intent(out) :: fluxes(mee_fluxes_nenergy)
         210             :     logical, intent(out) :: valid(mee_fluxes_nenergy)
         211             :
         212             :     integer :: i, ndx1, ndx2
         213             :     logical :: found
         214             :     real(r8) :: wght1,wght2
         215             :
         216           0 :     valid(:) = .not. mee_fluxes_fillin
         217           0 :     fluxes(:) = 0._r8
         218             :
         219           0 :     if (.not.mee_fluxes_active) return
         220             :
         221           0 :     found = .false.
         222             :
         223           0 :     findloop: do i = 1,nlshells-1
         224           0 :        if ( l_shell>=lshell(i) .and. l_shell<=lshell(i+1) ) then
         225           0 :           ndx1=i
         226           0 :           ndx2=i+1
         227           0 :           wght2 = (l_shell-lshell(ndx1))/(lshell(ndx2)-lshell(ndx1))
         228           0 :           wght1 = 1._r8 - wght2
         229           0 :           found = .true.
         230           0 :           exit findloop
         231             :        endif
         232             :     end do findloop
         233             :
         234           0 :     if (found) then
         235           0 :        if (mee_fluxes_fillin) then
         236           0 :           valid(:) = valflx(:,ndx1) .and. valflx(:,ndx2)
         237             :        end if
         238           0 :        where( valid(:) )
         239           0 :           fluxes(:) = wght1*influx(:,ndx1) + wght2*influx(:,ndx2)
         240             :        end where
         241             :     end if
         242             :
         243             :   end subroutine mee_fluxes_extract
         244             :
         245             :   !-----------------------------------------------------------------------------
         246             :   !-----------------------------------------------------------------------------
         247           0 :   subroutine read_fluxes()
         248             :
         249             :     ! local vars
         250             :     integer :: ierr, cnt(4), start(4)
         251             :
         252           0 :     cnt = (/1, mee_fluxes_nenergy, nlshells, 2/)
         253             :
         254             :     ! use the 50 percentile level data ( index 3 )
         255           0 :     start = (/3, 1, 1, time_coord%indxs(1)/)
         256             :
         257             :     ! float RBSP_flux_scaled(time, lshell, energy, percentiles) ;
         258           0 :     ierr = pio_get_var( file_id, flux_var_id, start, cnt, indata )
         259             :
         260           0 :   end subroutine read_fluxes
         261             :
         262             : end module mee_fluxes

Generated by: LCOV version 1.14