LCOV - code coverage report
Current view: | top level - hemco/HEMCO/src/Extensions - hcox_volcano_mod.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 0 | 394 | 0.0 %
Date: | 2024-12-17 17:57:11 |  | Functions: | 0 | 11 | 0.0 %

              Line data    Source code
           1             : !------------------------------------------------------------------------------
           2             : !                   Harmonized Emissions Component (HEMCO)                    !
           3             : !------------------------------------------------------------------------------
           4             : !BOP
           5             : !
           6             : ! !MODULE: hcox_volcano_mod.F90
           7             : !
           8             : ! !DESCRIPTION: Module HCOX\_Volcano\_Mod.F90 is a HEMCO extension to use
           9             : !  volcano emissions (such as AeroCom or OMI) from ascii tables. This module
          10             : !  reads the daily data tables and emits the emissions according to the
          11             : !  information in this file.
          12             : !\\
          13             : !\\
          14             : ! !INTERFACE:
          15             : !
          16             : MODULE HCOX_Volcano_Mod
          17             : !
          18             : ! !USES:
          19             : !
          20             :   USE HCO_Error_MOD
          21             :   USE HCO_Diagn_MOD
          22             :   USE HCOX_TOOLS_MOD
          23             :   USE HCOX_State_MOD, ONLY : Ext_State
          24             :   USE HCO_State_MOD,  ONLY : HCO_State
          25             :
          26             :   IMPLICIT NONE
          27             :   PRIVATE
          28             : !
          29             : ! !PUBLIC MEMBER FUNCTIONS:
          30             : !
          31             :   PUBLIC :: HCOX_Volcano_Run
          32             :   PUBLIC :: HCOX_Volcano_Init
          33             :   PUBLIC :: HCOX_Volcano_Final
          34             : !
          35             : ! !PRIVATE MEMBER FUNCTIONS:
          36             : !
          37             :   PRIVATE :: ReadVolcTable
          38             :   PRIVATE :: EmitVolc
          39             : !
          40             : ! !REMARKS:
          41             : ! Each volcano table (e.g. AeroCom or OMI) is expected to list the volcano
          42             : ! location, sulfur emissions (in kg S/s), and the volcano elevation as well
          43             : ! as the volcano plume column height. These entries need be separated by space
          44             : ! characters. For example:
          45             : !                                                                             .
          46             : ! ###  LAT (-90,90), LON (-180,180), SULFUR [kg S/s], ELEVATION [m], CLOUD_COLUMN_HEIGHT [m]
          47             : ! ### If elevation=cloud_column_height, emit in layer of elevation
          48             : ! ### else, emit in top 1/3 of cloud_column_height
          49             : ! volcano::
          50             : ! 50.170 6.850 3.587963e-03 600. 600.
          51             : ! ::
          52             : !                                                                             .
          53             : ! The sulfur read from table is emitted as the species defined in the
          54             : ! Volcano settings section. More than one species can be provided. Mass
          55             : ! sulfur is automatically converted to mass of emitted species (using the
          56             : ! emitted molecular weight and molecular ratio of the corresponding HEMCO
          57             : ! species). Additional scale factors can be defined in the settings section
          58             : ! by using the (optional) setting 'Scaling_<SpecName>'.
          59             : ! For example, to emit SO2 and BrO from volcanoes, with an additional scale
          60             : ! factor of 1e-4 kg BrO / kgS for BrO, use the following setting:
          61             : !
          62             : !117     Volcano                : on    SO2/BrO
          63             : !    --> Scaling_BrO            :       1.0e-4
          64             : !    --> Volcano_Source         :       AeroCom
          65             : !    --> Volcano_Table          :       $ROOT/VOLCANO/v2021-09/$YYYY/$MM/so2_volcanic_emissions_Carns.$YYYY$MM$DD.rc
          66             : !    --> Volcano_Climatology    :       $ROOT/VOLCANO/v2021-09/so2_volcanic_emissions_CARN_v202005.degassing_only.rc
          67             : !                                                                             .
          68             : ! This extension was originally added for usage within GEOS-5 and AeroCom
          69             : ! volcanic emissions, but has been modified to work with OMI-based volcanic
          70             : ! emissions from Ge et al. (2016).
          71             : !                                                                             .
          72             : ! When using this extension, you should turn off any other volcano emission
          73             : ! inventories!
          74             : !                                                                             .
          75             : !  References:
          76             : !  ============================================================================
          77             : !  (1 ) Ge, C., J. Wang, S. Carn, K. Yang, P. Ginoux, and N. Krotkov,
          78             : !       Satellite-based global volcanic SO2 emissions and sulfate direct
          79             : !       radiative forcing during 2005-2012, J. Geophys. Res. Atmos., 121(7),
          80             : !       3446-3464, doi:10.1002/2015JD023134, 2016.
          81             : !
          82             : ! !REVISION HISTORY:
          83             : !  04 Jun 2015 - C. Keller   - Initial version
          84             : !  See https://github.com/geoschem/hemco for complete history
          85             : !EOP
          86             : !------------------------------------------------------------------------------
          87             : !BOC
          88             : !
          89             : ! !MODULE VARIABLES:
          90             : !
          91             :   TYPE :: MyInst
          92             :    INTEGER                         :: Instance
          93             :    INTEGER                         :: ExtNr     = -1   ! Extension number
          94             :    INTEGER                         :: CatErupt  = -1   ! Category of eruptive emissions
          95             :    INTEGER                         :: CatDegas  = -1   ! Category of degassing emissions
          96             :    INTEGER                         :: nSpc      =  0   ! # of species
          97             :    INTEGER                         :: nVolc     =  0   ! # of volcanoes in buffer
          98             :    INTEGER,  ALLOCATABLE           :: SpcIDs(:)        ! HEMCO species IDs
          99             :    REAL(sp), ALLOCATABLE           :: SpcScl(:)        ! Species scale factors
         100             :    REAL(sp), ALLOCATABLE           :: VolcSlf(:)       ! Sulface emissions [kg S/s]
         101             :    REAL(sp), ALLOCATABLE           :: VolcElv(:)       ! Elevation [m]
         102             :    REAL(sp), ALLOCATABLE           :: VolcCld(:)       ! Cloud column height [m]
         103             :    INTEGER,  ALLOCATABLE           :: VolcIdx(:)       ! Lon grid index
         104             :    INTEGER,  ALLOCATABLE           :: VolcJdx(:)       ! Lat grid index
         105             :    INTEGER,  ALLOCATABLE           :: VolcBeg(:)       ! Begin time (optional)
         106             :    INTEGER,  ALLOCATABLE           :: VolcEnd(:)       ! End time   (optional)
         107             :    CHARACTER(LEN=255)              :: FileName         ! Volcano file name
         108             :    CHARACTER(LEN=255)              :: ClimFile         ! Climatology file name
         109             :    CHARACTER(LEN=255)              :: VolcSource       ! Volcano data source
         110             :    INTEGER                         :: YmdOnFile = -1   ! Date of file currently in record
         111             :    CHARACTER(LEN=61), ALLOCATABLE  :: SpcScalFldNme(:) ! Names of scale factor fields
         112             :    TYPE(MyInst), POINTER           :: NextInst => NULL()
         113             :   END TYPE MyInst
         114             :
         115             :   ! Pointer to instances
         116             :   TYPE(MyInst), POINTER            :: AllInst => NULL()
         117             :
         118             :   ! Volcano data is in kgS. Will be converted to kg emitted species.
         119             :   ! MW_S is the molecular weight of sulfur
         120             :   REAL(hp), PARAMETER             :: MW_S = 32.0_hp
         121             :
         122             : CONTAINS
         123             : !EOC
         124             : !------------------------------------------------------------------------------
         125             : !                   Harmonized Emissions Component (HEMCO)                    !
         126             : !------------------------------------------------------------------------------
         127             : !BOP
         128             : !
         129             : ! !IROUTINE: HCOX_Volcano_Run
         130             : !
         131             : ! !DESCRIPTION: Subroutine HCOX\_Volcano\_Run is the driver routine
         132             : ! for the customizable HEMCO extension.
         133             : !\\
         134             : !\\
         135             : ! !INTERFACE:
         136             : !
         137           0 :   SUBROUTINE HCOX_Volcano_Run( ExtState, HcoState, RC )
         138             : !
         139             : ! !USES:
         140             : !
         141             :     USE HCO_FluxArr_Mod,  ONLY : HCO_EmisAdd
         142             : !
         143             : ! !INPUT PARAMETERS:
         144             : !
         145             :     TYPE(Ext_State), POINTER       :: ExtState    ! Module options
         146             : !
         147             : ! !INPUT/OUTPUT PARAMETERS:
         148             : !
         149             :     TYPE(HCO_State), POINTER       :: HcoState    ! Hemco state
         150             :     INTEGER,         INTENT(INOUT) :: RC          ! Success or failure
         151             : !
         152             : ! !REVISION HISTORY:
         153             : !  04 Jun 2015 - C. Keller   - Initial version
         154             : !  See https://github.com/geoschem/hemco for complete history
         155             : !EOP
         156             : !------------------------------------------------------------------------------
         157             : !BOC
         158             : !
         159             : ! !LOCAL VARIABLES:
         160             : !
         161             :     ! Scalars
         162             :     INTEGER               :: N
         163             :     LOGICAL               :: ERR
         164             :
         165             :     ! Strings
         166             :     CHARACTER(LEN=255)    :: ErrMsg, ThisLoc, LOC
         167             :
         168             :     ! Arrays
         169           0 :     REAL(sp)              :: SO2degas(HcoState%NX,HcoState%NY,HcoState%NZ)
         170           0 :     REAL(sp)              :: SO2erupt(HcoState%NX,HcoState%NY,HcoState%NZ)
         171           0 :     REAL(sp)              :: iFlx    (HcoState%NX,HcoState%NY,HcoState%NZ)
         172             :
         173             :     ! Pointers
         174             :     TYPE(MyInst), POINTER :: Inst
         175             :
         176             :     !=================================================================
         177             :     ! HCOX_VOLCANO_RUN begins here!
         178             :     !=================================================================
         179           0 :     LOC = 'HCOX_VOLCANO_RUN (HCOX_VOLCANO_MOD.F90)'
         180             :
         181             :     ! Assume success
         182           0 :     RC = HCO_SUCCESS
         183             :
         184             :     ! Sanity check: return if extension not turned on
         185           0 :     IF ( ExtState%Volcano <= 0 ) RETURN
         186             :
         187             :     ! Enter
         188           0 :     CALL HCO_ENTER( HcoState%Config%Err, LOC, RC )
         189           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         190           0 :         CALL HCO_ERROR( 'ERROR 0', RC, THISLOC=LOC )
         191           0 :         RETURN
         192             :     ENDIF
         193             :
         194             :     ! Define strings for error messgaes
         195           0 :     ErrMsg = ''
         196             :     ThisLoc =  &
         197           0 :     ' -> in HCOX_Volcano_Run (in module HEMCO/Extensions/hcox_volcano_mod.F90)'
         198             :
         199             :     ! Get instance
         200           0 :     Inst => NULL()
         201           0 :     CALL InstGet( ExtState%Volcano, Inst, RC )
         202           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         203           0 :        WRITE( ErrMsg, * ) 'Cannot find Volcano instance Nr. ', ExtState%Volcano
         204           0 :        CALL HCO_Error( ErrMsg, RC, ThisLoc )
         205           0 :        RETURN
         206             :     ENDIF
         207             :
         208             :     !----------------------------------------------
         209             :     ! Read/update the volcano data
         210             :     ! (will be done only if this is a new day)
         211             :     !----------------------------------------------
         212           0 :     CALL ReadVolcTable( HcoState, ExtState, Inst, RC )
         213           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         214           0 :        ErrMsg = 'Error encountered in "ReadVolcTable"!'
         215           0 :        CALL HCO_Error( ErrMsg, RC, ThisLoc )
         216           0 :        RETURN
         217             :     ENDIF
         218             :
         219             :     !=======================================================================
         220             :     ! Compute volcano emissions for non dry-run simulations
         221             :     ! (Skip for GEOS-Chem dry-run or HEMCO-standalone dry-run)
         222             :     !=======================================================================
         223           0 :     IF ( .not. HcoState%Options%IsDryRun ) THEN
         224             :
         225             :        ! Emit volcanos into SO2degas and SO2erupt arrays [kg S/m2/s]
         226             :        CALL EmitVolc( HcoState, ExtState, Inst, &
         227           0 :                       SO2degas, SO2erupt, RC    )
         228           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         229           0 :           ErrMsg = 'Error encountered in "EmitVolc"!'
         230           0 :           CALL HCO_Error( ErrMsg, RC, ThisLoc )
         231           0 :           RETURN
         232             :        ENDIF
         233             :
         234             :        ! Add eruptive and degassing emissions to emission arrays & diagnostics
         235           0 :        DO N = 1, Inst%nSpc
         236             :
         237             :           !-------------------------------------------
         238             :           ! Add degassing emissions
         239             :           !-------------------------------------------
         240             :
         241             :           ! Convert from [kg S/m2/s] to [kg species/m2/s]
         242           0 :           iFlx = SO2degas * Inst%SpcScl(N)
         243             :
         244             :           ! Apply user-defined scaling (if any) for this species
         245             :           CALL HCOX_Scale( HcoState, iFlx, &
         246           0 :                            TRIM(Inst%SpcScalFldNme(N)), RC )
         247           0 :           IF ( RC /= HCO_SUCCESS ) THEN
         248           0 :              ErrMsg = 'Error encountered in "HCOX_Scale (degassing)"!'
         249           0 :              CALL HCO_Error( ErrMsg, RC, ThisLoc )
         250           0 :              RETURN
         251             :           ENDIF
         252             :
         253             :           ! Add degassing emissions into the HEMCO state
         254           0 :           CALL HCO_EmisAdd( HcoState, iFlx, Inst%SpcIDs(N), &
         255           0 :                             RC, ExtNr=Inst%ExtNr, Cat=Inst%CatDegas )
         256           0 :           IF ( RC /= HCO_SUCCESS ) THEN
         257           0 :              ErrMsg = 'Error encountered in "HCO_EmisAdd" (degassing)!'
         258           0 :              CALL HCO_Error( ErrMsg, RC, ThisLoc )
         259           0 :              RETURN
         260             :           ENDIF
         261             :
         262             :           !-------------------------------------------
         263             :           ! Add eruptive emissions
         264             :           !-------------------------------------------
         265             :
         266             :           ! Convert from [kg S/m2/s] to [kg species/m2/s]
         267           0 :           iFlx = SO2erupt * Inst%SpcScl(N)
         268             :
         269             :           ! Apply user-defined scaling (if any) for this species
         270             :           CALL HCOX_Scale( HcoState, iFlx, &
         271           0 :                            TRIM(Inst%SpcScalFldNme(N)), RC )
         272           0 :           IF ( RC /= HCO_SUCCESS ) THEN
         273           0 :              ErrMsg = 'Error encountered in "HCOX_Scale" (eruptive"!'
         274           0 :              CALL HCO_Error( ErrMsg, RC, ThisLoc )
         275           0 :              RETURN
         276             :           ENDIF
         277             :
         278             :           ! Add eruptive emissions to the HEMCO state
         279           0 :           CALL HCO_EmisAdd( HcoState, iFlx, Inst%SpcIDs(N), &
         280           0 :                             RC, ExtNr=Inst%ExtNr, Cat=Inst%CatErupt )
         281           0 :           IF ( RC /= HCO_SUCCESS ) THEN
         282           0 :              ErrMsg = 'Error encountered in "HCO_EmisAdd" (eruptive)!'
         283           0 :              CALL HCO_Error( ErrMsg, RC, ThisLoc )
         284           0 :              RETURN
         285             :           ENDIF
         286             :
         287             :        ENDDO !N
         288             :     ENDIF
         289             :
         290             :     !=======================================================================
         291             :     ! Exit
         292             :     !=======================================================================
         293             :
         294             :     ! Cleanup
         295           0 :     Inst => NULL()
         296             :
         297             :     ! Return w/ success
         298           0 :     CALL HCO_LEAVE( HcoState%Config%Err, RC )
         299             :
         300             :   END SUBROUTINE HCOX_Volcano_Run
         301             : !EOC
         302             : !------------------------------------------------------------------------------
         303             : !                   Harmonized Emissions Component (HEMCO)                    !
         304             : !------------------------------------------------------------------------------
         305             : !BOP
         306             : !
         307             : ! !IROUTINE: HCOX_Volcano_Init
         308             : !
         309             : ! !DESCRIPTION: Subroutine HCOX\_Volcano\_Init initializes the HEMCO
         310             : ! CUSTOM extension.
         311             : !\\
         312             : !\\
         313             : ! !INTERFACE:
         314             : !
         315           0 :   SUBROUTINE HCOX_Volcano_Init( HcoState, ExtName,ExtState, RC )
         316             : !
         317             : ! !USES:
         318             : !
         319             :     USE HCO_ExtList_Mod,    ONLY : GetExtNr
         320             :     USE HCO_ExtList_Mod,    ONLY : GetExtOpt
         321             :     USE HCO_ExtList_Mod,    ONLY : GetExtSpcVal
         322             :     USE HCO_STATE_MOD,      ONLY : HCO_GetExtHcoID
         323             : !
         324             : ! !INPUT PARAMETERS:
         325             : !
         326             :     CHARACTER(LEN=*), INTENT(IN   ) :: ExtName    ! Extension name
         327             :     TYPE(Ext_State),  POINTER       :: ExtState   ! Module options
         328             : !
         329             : ! !INPUT/OUTPUT PARAMETERS:
         330             : !
         331             :     TYPE(HCO_State),  POINTER       :: HcoState   ! Hemco state
         332             :     INTEGER,          INTENT(INOUT) :: RC
         333             :
         334             : ! !REVISION HISTORY:
         335             : !  04 Jun 2015 - C. Keller   - Initial version
         336             : !  See https://github.com/geoschem/hemco for complete history
         337             : !EOP
         338             : !------------------------------------------------------------------------------
         339             : !BOC
         340             : !
         341             : ! !LOCAL VARIABLES:
         342             : !
         343             :     TYPE(MyInst), POINTER          :: Inst
         344             :     REAL(sp)                       :: ValSp
         345             :     INTEGER                        :: ExtNr, N, Dum
         346             :     LOGICAL                        :: FOUND
         347           0 :     CHARACTER(LEN=31), ALLOCATABLE :: SpcNames(:)
         348             :     CHARACTER(LEN=255)             :: MSG, Str, LOC
         349             :
         350             :     !=================================================================
         351             :     ! HCOX_VOLCANO_INIT begins here!
         352             :     !=================================================================
         353           0 :     LOC = 'HCOX_VOLCANO_INIT (HCOX_VOLCANO_MOD.F90)'
         354             :
         355             :     ! Extension Nr.
         356           0 :     ExtNr = GetExtNr( HcoState%Config%ExtList, TRIM(ExtName) )
         357           0 :     IF ( ExtNr <= 0 ) RETURN
         358             :
         359             :     ! Enter
         360           0 :     CALL HCO_ENTER( HcoState%Config%Err, LOC, RC )
         361           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         362           0 :         CALL HCO_ERROR( 'ERROR 1', RC, THISLOC=LOC )
         363           0 :         RETURN
         364             :     ENDIF
         365             :
         366             :     ! Create instance for this simulation
         367           0 :     Inst => NULL()
         368           0 :     CALL InstCreate( ExtNr, ExtState%Volcano, Inst, RC )
         369           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         370           0 :        CALL HCO_Error( 'Cannot create Volcano instance', RC )
         371           0 :        RETURN
         372             :     ENDIF
         373             :
         374             :     ! Write the name of the extension regardless of the verbose settings
         375           0 :     IF ( HcoState%amIRoot ) THEN
         376           0 :        msg = 'Using HEMCO extension: Volcano (volcanic SO2 emissions)'
         377           0 :        IF ( HCO_IsVerb( HcoState%Config%Err ) ) THEN
         378           0 :           CALL HCO_Msg( HcoState%Config%Err, msg, sep1='-' ) ! with separator
         379             :        ELSE
         380           0 :           CALL HCO_Msg( msg, verb=.TRUE.                   ) ! w/o separator
         381             :        ENDIF
         382             :     ENDIF
         383             :
         384             :     ! Get species IDs.
         385             :     CALL HCO_GetExtHcoID( HcoState, ExtNr,     Inst%SpcIDs,                  &
         386           0 :                           SpcNames, Inst%nSpc, RC                           )
         387           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         388           0 :         CALL HCO_ERROR( 'ERROR 2', RC, THISLOC=LOC )
         389           0 :         RETURN
         390             :     ENDIF
         391             :
         392             :     ! There must be at least one species
         393           0 :     IF ( Inst%nSpc == 0 ) THEN
         394             :        CALL HCO_Error(                                  &
         395           0 :                       'No Volcano species specified', RC                    )
         396           0 :        RETURN
         397             :     ENDIF
         398             :
         399             :     ! Determine scale factor to be applied to each species. This is 1.00
         400             :     ! by default, but can be set in the HEMCO configuration file via setting
         401             :     ! Scaling_<SpcName>.
         402             :     CALL GetExtSpcVal( HcoState%Config, ExtNr,  Inst%nSpc,   SpcNames,       &
         403           0 :                        'Scaling',       1.0_sp, Inst%SpcScl, RC             )
         404           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         405           0 :         CALL HCO_ERROR( 'ERROR 3', RC, THISLOC=LOC )
         406           0 :         RETURN
         407             :     ENDIF
         408             :
         409             :     ! Get species mask fields
         410             :     CALL GetExtSpcVal( HcoState%Config,    ExtNr,        Inst%nSpc,          &
         411             :                        SpcNames,           'ScaleField', HCOX_NOSCALE,       &
         412           0 :                        Inst%SpcScalFldNme, RC                               )
         413           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         414           0 :         CALL HCO_ERROR( 'ERROR 4', RC, THISLOC=LOC )
         415           0 :         RETURN
         416             :     ENDIF
         417             :
         418             :     ! Add conversion factor from kg S to kg species
         419           0 :     DO N = 1, Inst%nSpc
         420           0 :        Inst%SpcScl(N) = Inst%SpcScl(N) * HcoState%Spc(Inst%SpcIDs(N))%MW_g &
         421           0 :                         / MW_S
         422             :     ENDDO
         423             :
         424             :     ! Get location of volcano table. This must be provided.
         425             :     CALL GetExtOpt( HcoState%Config, ExtNr, 'Volcano_Table',                 &
         426           0 :                     OptValChar=Inst%FileName, FOUND=FOUND, RC=RC            )
         427             :
         428           0 :     IF ( RC /= HCO_SUCCESS .OR. .NOT. FOUND ) THEN
         429             :        MSG = 'Cannot read Volcano table file name. Please provide '       // &
         430             :              'the Volcano table as a setting to the Volcano extension. '  // &
         431           0 :              'The name of this setting must be `Volcano_Table`.'
         432           0 :        CALL HCO_Error( MSG, RC )
         433           0 :        RETURN
         434             :     ENDIF
         435             :
         436             :     ! Get location of volcano climatology table. This must be provided.
         437             :     CALL GetExtOpt( HcoState%Config, ExtNr, 'Volcano_Climatology',           &
         438           0 :                     OptValChar=Inst%ClimFile, FOUND=FOUND, RC=RC            )
         439             :
         440           0 :     IF ( RC /= HCO_SUCCESS .OR. .NOT. FOUND ) THEN
         441             :        MSG = 'Cannot read Volcano climatology file name. Please provide ' // &
         442             :              'the Volcano climatology as a setting to the Volcano extension. ' // &
         443           0 :              'The name of this setting must be `Volcano_Climatology`.'
         444           0 :        CALL HCO_Error( HcoState%Config%Err, MSG, RC )
         445           0 :        RETURN
         446             :     ENDIF
         447             :
         448             :     ! See if emissions data source is given
         449             :     ! As of v11-02f, options are AeroCom or OMI
         450           0 :     Inst%VolcSource = 'AeroCom'
         451             :     CALL GetExtOpt( HcoState%Config, ExtNr,       'Volcano_Source',          &
         452           0 :                     OptValChar=Str,  FOUND=FOUND, RC=RC                     )
         453           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         454           0 :         CALL HCO_ERROR( 'ERROR 5', RC, THISLOC=LOC )
         455           0 :         RETURN
         456             :     ENDIF
         457           0 :     IF ( FOUND ) Inst%VolcSource = Str
         458             :
         459             :     ! See if eruptive and degassing hierarchies are given
         460           0 :     Inst%CatErupt = 51
         461           0 :     Inst%CatDegas = 52
         462             :     CALL GetExtOpt( HcoState%Config, ExtNr,       'Cat_Degassing',           &
         463           0 :                     OptValInt=Dum,   FOUND=FOUND, RC=RC                     )
         464           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         465           0 :         CALL HCO_ERROR( 'ERROR 6', RC, THISLOC=LOC )
         466           0 :         RETURN
         467             :     ENDIF
         468           0 :     IF ( FOUND ) Inst%CatDegas = Dum
         469             :     CALL GetExtOpt( HcoState%Config, ExtNr,       'Cat_Eruptive',            &
         470           0 :                     OptValInt=Dum,   FOUND=FOUND, RC=RC                    )
         471           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         472           0 :         CALL HCO_ERROR( 'ERROR 7', RC, THISLOC=LOC )
         473           0 :         RETURN
         474             :     ENDIF
         475           0 :     IF ( FOUND ) Inst%CatErupt = Dum
         476             :
         477             :     ! Verbose mode
         478           0 :     IF ( HcoState%amIRoot ) THEN
         479           0 :        MSG = ' - use the following species (Name, HcoID, Scaling relative to kgS):'
         480           0 :        CALL HCO_MSG( HcoState%Config%Err, MSG)
         481           0 :        DO N = 1, Inst%nSpc
         482           0 :           WRITE(MSG,*) TRIM(SpcNames(N)), ', ', Inst%SpcIDs(N), ', ', Inst%SpcScl(N)
         483           0 :           CALL HCO_MSG( HcoState%Config%Err, MSG)
         484           0 :           WRITE(MSG,*) 'Apply scale field: ', TRIM(Inst%SpcScalFldNme(N))
         485           0 :           CALL HCO_MSG( HcoState%Config%Err, MSG)
         486             :        ENDDO
         487           0 :        WRITE(MSG,*) ' - Emissions data source is ', TRIM(Inst%VolcSource)
         488           0 :        CALL HCO_MSG( HcoState%Config%Err,  MSG )
         489           0 :        WRITE(MSG,*) ' - Emit eruptive emissions as category ', Inst%CatErupt
         490           0 :        CALL HCO_MSG( HcoState%Config%Err,  MSG )
         491           0 :        WRITE(MSG,*) ' - Emit degassing emissions as category ', Inst%CatDegas
         492           0 :        CALL HCO_MSG( HcoState%Config%Err,  MSG )
         493             :     ENDIF
         494             :
         495             :     ! Cleanup
         496           0 :     Inst => NULL()
         497           0 :     IF ( ALLOCATED(SpcNames) ) DEALLOCATE(SpcNames)
         498             :
         499           0 :     CALL HCO_Leave( HcoState%Config%Err, RC )
         500             :
         501           0 :   END SUBROUTINE HCOX_Volcano_Init
         502             : !EOC
         503             : !------------------------------------------------------------------------------
         504             : !                   Harmonized Emissions Component (HEMCO)                    !
         505             : !------------------------------------------------------------------------------
         506             : !BOP
         507             : !
         508             : ! !IROUTINE: HCOX_Volcano_Final
         509             : !
         510             : ! !DESCRIPTION: Subroutine HCOX\_AeroCom\_Final finalizes the HEMCO
         511             : !  AeroCom extension.
         512             : !\\
         513             : !\\
         514             : ! !INTERFACE:
         515             : !
         516           0 :   SUBROUTINE HCOX_Volcano_Final( ExtState )
         517             : !
         518             : ! !INPUT PARAMETERS:
         519             : !
         520             :     TYPE(Ext_State), POINTER :: ExtState   ! Module options
         521             : !
         522             : ! !REVISION HISTORY:
         523             : !  04 Jun 2015 - C. Keller   - Initial version
         524             : !  See https://github.com/geoschem/hemco for complete history
         525             : !EOP
         526             : !------------------------------------------------------------------------------
         527             : !BOC
         528             :     !=================================================================
         529             :     ! HCOX_VOLCANO_FINAL begins here!
         530             :     !=================================================================
         531           0 :     CALL InstRemove( ExtState%Volcano )
         532             :
         533           0 :   END SUBROUTINE HCOX_Volcano_Final
         534             : !EOC
         535             : !------------------------------------------------------------------------------
         536             : !                   Harmonized Emissions Component (HEMCO)                    !
         537             : !------------------------------------------------------------------------------
         538             : !BOP
         539             : !
         540             : ! !IROUTINE: ReadVolcTable
         541             : !
         542             : ! !DESCRIPTION: Subroutine ReadVolcTable reads the AeroCom volcano table of the
         543             : !  current day.
         544             : !\\
         545             : !\\
         546             : ! !INTERFACE:
         547             : !
         548           0 :   SUBROUTINE ReadVolcTable( HcoState, ExtState, Inst, RC )
         549             : !
         550             : ! !USES:
         551             : !
         552             :     USE HCO_CharTools_Mod
         553             :     USE HCO_inquireMod,     ONLY : findfreeLun
         554             :     USE HCO_CLOCK_MOD,      ONLY : HcoClock_NewDay
         555             :     USE HCO_CLOCK_MOD,      ONLY : HcoClock_Get
         556             :     USE HCO_GeoTools_MOD,   ONLY : HCO_GetHorzIJIndex
         557             :     USE HCO_EXTLIST_MOD,    ONLY : HCO_GetOpt
         558             : !
         559             : ! !INPUT PARAMETERS:
         560             : !
         561             :     TYPE(Ext_State),  POINTER       :: ExtState   ! Module options
         562             : !
         563             : ! !INPUT/OUTPUT PARAMETERS:
         564             : !
         565             :     TYPE(HCO_State),  POINTER       :: HcoState   ! Hemco state
         566             :     TYPE(MyInst),     POINTER       :: Inst
         567             :     INTEGER,          INTENT(INOUT) :: RC
         568             :
         569             : ! !REVISION HISTORY:
         570             : !  04 Jun 2015 - C. Keller   - Initial version
         571             : !  See https://github.com/geoschem/hemco for complete history
         572             : !EOP
         573             : !------------------------------------------------------------------------------
         574             : !BOC
         575             : !
         576             : ! !LOCAL VARIABLES:
         577             : !
         578             :     INTEGER               :: YYYY, MM, DD
         579             :     INTEGER               :: ThisYMD
         580             :     INTEGER               :: N, LUN, IOS, AS
         581             :     INTEGER               :: nVolc, nCol
         582             :     REAL(sp)              :: Dum(10)
         583           0 :     REAL(hp), ALLOCATABLE :: VolcLon(:)      ! Volcano longitude [deg E]
         584           0 :     REAL(hp), ALLOCATABLE :: VolcLat(:)      ! Volcano latitude  [deg N]
         585             :     LOGICAL               :: FileExists, EOF
         586             :     CHARACTER(LEN=255)    :: ThisFile, ThisLine
         587             :     CHARACTER(LEN=255)    :: MSG,      FileMsg
         588             :     CHARACTER(LEN=255)    :: LOC = 'ReadVolcTable (hcox_volcano_mod.F90)'
         589             :
         590             :     !=================================================================
         591             :     ! ReadVolcTable begins here!
         592             :     !=================================================================
         593             :
         594             :     ! Get current year, month, day
         595           0 :     CALL HcoClock_Get ( HcoState%Clock, cYYYY=YYYY, cMM=MM, cDD=DD, RC=RC )
         596           0 :     IF ( RC /= HCO_SUCCESS ) RETURN
         597             : #if defined( MODEL_GEOS )
         598             :     ! Error trap: skip leap days
         599             :     IF ( MM == 2 .AND. DD > 28 ) DD = 28
         600             : #endif
         601             :
         602             :     ! Compare current day against day on file
         603           0 :     ThisYMD  = YYYY*10000 + MM*100+ DD
         604             :
         605             :     ! Do only if it's a different day
         606           0 :     IF ( ThisYMD /= Inst%YmdOnFile ) THEN
         607             :
         608             :        ! Get file name
         609           0 :        ThisFile = Inst%FileName
         610           0 :        CALL HCO_CharParse( HcoState%Config, ThisFile, YYYY, MM, DD, 0, 0, RC )
         611           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         612           0 :            CALL HCO_ERROR( 'ERROR 9', RC, THISLOC=LOC )
         613           0 :            RETURN
         614             :        ENDIF
         615             :
         616             :        !--------------------------------------------------------------------
         617             :        ! In dry-run mode, print file path to dryrun log and exit.
         618             :        ! Otherwise, print file path to the HEMCO log file and continue.
         619             :        !--------------------------------------------------------------------
         620             :
         621             :        ! Test if the file exists
         622           0 :        INQUIRE( FILE=TRIM( ThisFile ), EXIST=FileExists )
         623             :
         624             :        ! Create a display string based on whether or not the file is found
         625           0 :        IF ( FileExists ) THEN
         626           0 :           FileMsg = 'HEMCO (VOLCANO): Opening'
         627             :        ELSE
         628           0 :           FileMsg = 'HEMCO (VOLCANO): REQUIRED FILE NOT FOUND'
         629             :        ENDIF
         630             :
         631             :        ! Write file status to stdout and the HEMCO log
         632           0 :        IF ( Hcostate%amIRoot ) THEN
         633           0 :           WRITE( 6,   300 ) TRIM( FileMsg ), TRIM( ThisFile )
         634           0 :           WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( ThisFile )
         635           0 :           CALL HCO_MSG( HcoState%Config%Err, MSG )
         636             :  300      FORMAT( a, ' ', a )
         637             :        ENDIF
         638             :
         639           0 :        IF ( .not. FileExists ) THEN
         640             :
         641             :           ! Attempt to use climatology file instead
         642           0 :           ThisFile = Inst%ClimFile
         643             :           CALL HCO_CharParse( HcoState%Config, ThisFile, &
         644           0 :                               YYYY, MM, DD, 0, 0, RC )
         645           0 :           IF ( RC /= HCO_SUCCESS ) THEN
         646           0 :               CALL HCO_ERROR( 'ERROR 10', RC, THISLOC=LOC )
         647           0 :               RETURN
         648             :           ENDIF
         649             :
         650             :           ! Test if the file exists
         651           0 :           INQUIRE( FILE=TRIM( ThisFile ), EXIST=FileExists )
         652             :
         653             :           ! Write message to stdout and HEMCO log
         654           0 :           IF ( Hcostate%amIRoot ) THEN
         655           0 :              MSG = 'Attempting to read volcano climatology file'
         656           0 :              WRITE( 6,   300 ) TRIM( MSG )
         657           0 :              CALL HCO_MSG( HcoState%Config%Err, MSG )
         658             :           ENDIF
         659             :
         660             :           ! Create a display string based on whether or not the file is found
         661           0 :           IF ( FileExists ) THEN
         662           0 :              FileMsg = 'HEMCO (VOLCANO): Opening'
         663             :           ELSE
         664           0 :              FileMsg = 'HEMCO (VOLCANO): CLIMATOLOGY FILE NOT FOUND'
         665             :           ENDIF
         666             :
         667             :           ! Write file status to stdout and the HEMCO log
         668           0 :           IF ( Hcostate%amIRoot ) THEN
         669           0 :              WRITE( 6,   300 ) TRIM( FileMsg ), TRIM( ThisFile )
         670           0 :              WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( ThisFile )
         671           0 :              CALL HCO_MSG( HcoState%Config%Err, MSG )
         672             :           ENDIF
         673             :
         674             :        ENDIF
         675             :
         676             :        ! For dry-run simulations, return to calling program.
         677             :        ! For regular simulations, throw an error if we can't find the file.
         678           0 :        IF ( HcoState%Options%IsDryRun ) THEN
         679             :           RETURN
         680             :        ELSE
         681           0 :           IF ( .not. FileExists ) THEN
         682           0 :              WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( ThisFile )
         683           0 :              CALL HCO_ERROR( MSG, RC )
         684           0 :              RETURN
         685             :           ENDIF
         686             :        ENDIF
         687             :
         688             :        !--------------------------------------------------------------------
         689             :        ! Read data from files
         690             :        !--------------------------------------------------------------------
         691             :
         692             :        ! Open file
         693           0 :        LUN = findFreeLun()
         694           0 :        OPEN ( LUN, FILE=TRIM(ThisFile), STATUS='OLD', IOSTAT=IOS )
         695           0 :        IF ( IOS /= 0 ) THEN
         696           0 :           MSG = 'Error reading ' // TRIM(ThisFile)
         697           0 :           CALL HCO_ERROR(  MSG, RC, THISLOC=LOC )
         698           0 :           RETURN
         699             :        ENDIF
         700             :
         701             :        ! Get number of volcano records
         702           0 :        nVolc = 0
         703             :        DO
         704           0 :           CALL GetNextLine( LUN, ThisLine, EOF, RC )
         705           0 :           IF ( RC /= HCO_SUCCESS ) THEN
         706           0 :               CALL HCO_ERROR( 'ERROR 11', RC, THISLOC=LOC )
         707           0 :               RETURN
         708             :           ENDIF
         709           0 :           IF ( EOF ) EXIT
         710             :
         711             :           ! Skip any entries that contain '::'
         712           0 :           IF ( INDEX( TRIM(ThisLine), '::') > 0 ) CYCLE
         713             :
         714             :           ! If we make it to here, this is a valid entry
         715           0 :           nVolc = nVolc + 1
         716             :        ENDDO
         717             :
         718             :        ! Verbose
         719           0 :        IF ( HCO_IsVerb( HcoState%Config%Err ) ) THEN
         720           0 :           WRITE(MSG,*) 'Number of volcanoes: ', nVolc
         721           0 :           CALL HCO_MSG( HcoState%Config%Err, MSG)
         722             :        ENDIF
         723             :
         724             :        ! Allocate arrays
         725           0 :        IF ( nVolc > 0 ) THEN
         726             :           ! Eventually deallocate previously allocated data
         727           0 :           IF ( ALLOCATED(Inst%VolcSlf) ) DEALLOCATE(Inst%VolcSlf)
         728           0 :           IF ( ALLOCATED(Inst%VolcElv) ) DEALLOCATE(Inst%VolcElv)
         729           0 :           IF ( ALLOCATED(Inst%VolcCld) ) DEALLOCATE(Inst%VolcCld)
         730           0 :           IF ( ALLOCATED(Inst%VolcIdx) ) DEALLOCATE(Inst%VolcIdx)
         731           0 :           IF ( ALLOCATED(Inst%VolcJdx) ) DEALLOCATE(Inst%VolcJdx)
         732           0 :           IF ( ALLOCATED(Inst%VolcBeg) ) DEALLOCATE(Inst%VolcBeg)
         733           0 :           IF ( ALLOCATED(Inst%VolcEnd) ) DEALLOCATE(Inst%VolcEnd)
         734             :
         735             :           ALLOCATE(     VolcLon(nVolc), &
         736             :                         VolcLat(nVolc), &
         737             :                    Inst%VolcSlf(nVolc), &
         738             :                    Inst%VolcElv(nVolc), &
         739             :                    Inst%VolcCld(nVolc), &
         740             :                    Inst%VolcIdx(nVolc), &
         741             :                    Inst%VolcJdx(nVolc), &
         742             :                    Inst%VolcBeg(nVolc), &
         743             :                    Inst%VolcEnd(nVolc), &
         744           0 :                    STAT=AS )
         745           0 :           IF ( AS /= 0 ) THEN
         746             :              CALL HCO_ERROR ( &
         747           0 :                               'Volc allocation error', RC, THISLOC=LOC )
         748           0 :              RETURN
         749             :           ENDIF
         750           0 :                VolcLon = 0.0_hp
         751           0 :                VolcLat = 0.0_hp
         752           0 :           Inst%VolcSlf = 0.0_sp
         753           0 :           Inst%VolcElv = 0.0_sp
         754           0 :           Inst%VolcCld = 0.0_sp
         755           0 :           Inst%VolcBeg = 0
         756           0 :           Inst%VolcEnd = 0
         757             :
         758             :        ELSE
         759           0 :           WRITE(MSG,*) 'No volcano data found for year/mm/dd: ', YYYY, MM, DD
         760           0 :           CALL HCO_WARNING( HcoState%Config%Err, MSG, RC, THISLOC=LOC )
         761             :        ENDIF
         762             :
         763             :        ! Now read records
         764           0 :        IF ( nVolc > 0 ) THEN
         765           0 :           REWIND( LUN )
         766             :
         767           0 :           N = 0
         768             :           DO
         769           0 :              CALL GetNextLine( LUN, ThisLine, EOF, RC )
         770           0 :              IF ( RC /= HCO_SUCCESS ) THEN
         771           0 :                  CALL HCO_ERROR( 'ERROR 12', RC, THISLOC=LOC )
         772           0 :                  RETURN
         773             :              ENDIF
         774           0 :              IF ( EOF ) EXIT
         775             :
         776             :              ! Skip any entries that contain '::'
         777           0 :              IF ( INDEX( TRIM(ThisLine), '::') > 0 ) CYCLE
         778             :
         779             :              ! Write this data into the following vector element
         780           0 :              N = N + 1
         781           0 :              IF ( N > nVolc ) THEN
         782           0 :                 WRITE(MSG,*) 'N exceeds nVolc: ', N, nVolc, &
         783           0 :                              ' - This error occurred when reading ', &
         784           0 :                              TRIM(ThisFile), '. This line: ', TRIM(ThisLine)
         785           0 :                 CALL HCO_ERROR ( MSG, RC, THISLOC = LOC )
         786           0 :                 RETURN
         787             :              ENDIF
         788             :
         789             :              CALL HCO_CharSplit( TRIM(ThisLine), ' ', &
         790             :                                  HCO_GetOpt(HcoState%Config%ExtList,'Wildcard'), &
         791           0 :                                  Dum, nCol, RC )
         792           0 :              IF ( RC /= HCO_SUCCESS ) THEN
         793           0 :                  CALL HCO_ERROR( 'ERROR 13', RC, THISLOC=LOC )
         794           0 :                  RETURN
         795             :              ENDIF
         796             :
         797             :              ! Allow for 5 or 7 values
         798           0 :              IF ( nCol /= 5 .and. nCol /= 7 ) THEN
         799           0 :                 WRITE(MSG,*) 'Cannot parse line ', TRIM(ThisLine), &
         800           0 :                              'Expected five or seven entries, separated by ', &
         801           0 :                              'space character, instead found ', nCol
         802           0 :                 CALL HCO_ERROR ( MSG, RC, THISLOC = LOC )
         803           0 :                 RETURN
         804             :              ENDIF
         805             :
         806             :              ! Now pass to vectors
         807           0 :                   VolcLat(N) = Dum(1)
         808           0 :                   VolcLon(N) = Dum(2)
         809           0 :              Inst%VolcSlf(N) = Dum(3)
         810           0 :              Inst%VolcElv(N) = Dum(4)
         811           0 :              Inst%VolcCld(N) = Dum(5)
         812             :
         813             :              ! Some lines also include start and end time
         814           0 :              IF ( nCol == 7 ) THEN
         815           0 :                 Inst%VolcBeg(N) = Dum(6)
         816           0 :                 Inst%VolcEnd(N) = DUM(7)
         817             :              ENDIF
         818             :           ENDDO
         819             :
         820             :           ! At this point, we should have read exactly nVolc entries!
         821           0 :           IF ( N /= nVolc ) THEN
         822           0 :              WRITE(MSG,*) 'N /= nVolc: ', N, nVolc, &
         823           0 :                           ' - This error occurred when reading ', TRIM(ThisFile)
         824           0 :              CALL HCO_ERROR ( MSG, RC, THISLOC = LOC )
         825           0 :              RETURN
         826             :           ENDIF
         827             :
         828             :        ENDIF
         829             :
         830             :        ! All done
         831           0 :        CLOSE ( LUN )
         832             :
         833             :        ! Get grid box indeces for each location
         834           0 :        IF ( nVolc > 0 ) THEN
         835             :           CALL HCO_GetHorzIJIndex( HcoState, nVolc, VolcLon, &
         836           0 :                                    VolcLat, Inst%VolcIdx, Inst%VolcJdx, RC )
         837           0 :           IF ( RC /= HCO_SUCCESS ) THEN
         838           0 :               CALL HCO_ERROR( 'ERROR 14', RC, THISLOC=LOC )
         839           0 :               RETURN
         840             :           ENDIF
         841             :        ENDIF
         842             :
         843             :        ! Save # of volcanoes in archive
         844           0 :        Inst%nVolc = nVolc
         845             :
         846             :        ! Update date for file on record
         847           0 :        Inst%YmdOnFile = ThisYMD
         848             :
         849             :     ENDIF ! new day
         850             :
         851             :     ! Cleanup
         852           0 :     IF ( ALLOCATED(VolcLon) ) DEALLOCATE(VolcLon)
         853           0 :     IF ( ALLOCATED(VolcLat) ) DEALLOCATE(VolcLat)
         854             :
         855             :     ! Return w/ success
         856           0 :     RC = HCO_SUCCESS
         857             :
         858           0 :   END SUBROUTINE ReadVolcTable
         859             : !EOC
         860             : !------------------------------------------------------------------------------
         861             : !                   Harmonized Emissions Component (HEMCO)                    !
         862             : !------------------------------------------------------------------------------
         863             : !BOP
         864             : !
         865             : ! !IROUTINE: EmitVolc
         866             : !
         867             : ! !DESCRIPTION: Subroutine EmitVolc reads the AeroCom volcano table of the
         868             : !  current day.
         869             : !\\
         870             : !\\
         871             : ! !INTERFACE:
         872             : !
         873           0 :   SUBROUTINE EmitVolc( HcoState, ExtState, Inst, SO2d, SO2e, RC )
         874             : !
         875             : ! !USES:
         876             : !
         877             :     USE HCO_CLOCK_MOD,      ONLY : HcoClock_Get
         878             : !
         879             : ! !INPUT PARAMETERS:
         880             : !
         881             :     TYPE(Ext_State),  POINTER       :: ExtState   ! Module options
         882             : !
         883             : ! !INPUT/OUTPUT PARAMETERS:
         884             : !
         885             :     TYPE(HCO_State),  POINTER       :: HcoState   ! Hemco state
         886             :     TYPE(MyInst),     POINTER       :: Inst
         887             :     INTEGER,          INTENT(INOUT) :: RC
         888             : !
         889             : ! !OUTPUT PARAMETERS:
         890             : !
         891             :     REAL(sp),         INTENT(  OUT) :: SO2e(HcoState%NX,HcoState%NY,HcoState%NZ)
         892             :     REAL(sp),         INTENT(  OUT) :: SO2d(HcoState%NX,HcoState%NY,HcoState%NZ)
         893             : !
         894             : ! !REVISION HISTORY:
         895             : !  04 Jun 2015 - C. Keller   - Initial version
         896             : !  See https://github.com/geoschem/hemco for complete history
         897             : !EOP
         898             : !------------------------------------------------------------------------------
         899             : !BOC
         900             : !
         901             : ! !LOCAL VARIABLES:
         902             : !
         903             :     INTEGER            :: I, J, L, N, HH, MN, hhmmss
         904             :     LOGICAL            :: Erupt
         905             :     REAL(sp)           :: nSO2, zTop, zBot, PlumeHgt
         906             :     REAL(sp)           :: z1,   z2
         907             :     REAL(sp)           :: tmp1, tmp2, Frac
         908             :     REAL(sp)           :: totE, totD, volcE, volcD
         909             :     CHARACTER(LEN=255) :: MSG
         910             :     CHARACTER(LEN=255) :: LOC = 'EmitVolc (hcox_volcano_mod.F90)'
         911             :
         912             :     !=================================================================
         913             :     ! EmitVolc begins here!
         914             :     !=================================================================
         915             :
         916             :     ! Init
         917           0 :     SO2e = 0.0_sp
         918           0 :     SO2d = 0.0_sp
         919           0 :     totE = 0.0_sp
         920           0 :     totD = 0.0_sp
         921             :
         922             :     ! Make sure all required grid quantities are defined
         923           0 :     IF ( .NOT. ASSOCIATED(HcoState%Grid%AREA_M2%Val) ) THEN
         924             :        CALL HCO_ERROR ( &
         925           0 :                        'Grid box areas not defined', RC, THISLOC=LOC )
         926           0 :        RETURN
         927             :     ENDIF
         928           0 :     IF ( .NOT. ASSOCIATED(HcoState%Grid%ZSFC%Val) ) THEN
         929             :        CALL HCO_ERROR ( &
         930           0 :                        'Surface heights not defined', RC, THISLOC=LOC )
         931           0 :        RETURN
         932             :     ENDIF
         933           0 :     IF ( .NOT. ASSOCIATED(HcoState%Grid%BXHEIGHT_M%Val) ) THEN
         934             :        CALL HCO_ERROR ( &
         935           0 :                        'Grid box heights not defined', RC, THISLOC=LOC )
         936           0 :        RETURN
         937             :     ENDIF
         938             :
         939             :     ! Get current hour, minute and save as hhmmss
         940           0 :     CALL HcoClock_Get ( HcoState%Clock, cH=HH, cM=MN, RC=RC )
         941           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         942           0 :         CALL HCO_ERROR( 'ERROR 15', RC, THISLOC=LOC )
         943           0 :         RETURN
         944             :     ENDIF
         945           0 :     hhmmss = HH*10000 + MN*100
         946             :
         947             :     ! Do for every volcano
         948           0 :     IF ( Inst%nVolc > 0 ) THEN
         949           0 :        DO N = 1, Inst%nVolc
         950             :
         951             :           ! Grid box index for this volcano
         952           0 :           I = Inst%VolcIdx(N)
         953           0 :           J = Inst%VolcJdx(N)
         954             :
         955             :           ! Skip if outside of domain
         956           0 :           IF( I < 1 .OR. J < 1 ) CYCLE
         957             :
         958             :           ! Check if beginning and end time are specified
         959             :           ! Do not include emissions for this volcano outside of
         960             :           ! start and end times (mps, 6/20/19)
         961           0 :           IF ( Inst%VolcBeg(N) > 0 .or. Inst%VolcEnd(N) > 0 ) THEN
         962           0 :              IF ( hhmmss <  Inst%VolcBeg(N) ) CYCLE
         963           0 :              IF ( hhmmss >= Inst%VolcEnd(N) ) CYCLE
         964             :           ENDIF
         965             :
         966             :           ! total emissions of this volcano
         967           0 :           volcE = 0.0_sp
         968           0 :           volcD = 0.0_sp
         969             :
         970           0 :           z1 = HcoState%Grid%ZSFC%Val(I,J)
         971             :
         972             :           ! Get total emitted kgS/m2/s. Data in table is in kgS/s.
         973           0 :           nSo2 = Inst%VolcSlf(N) / HcoState%Grid%AREA_M2%Val(I,J)
         974             :
         975             :           ! Elevation of volcano base and volcano cloud top height [m]
         976             :           ! Make sure that the bottom / top are at least at surface level
         977           0 :           zBot = MAX(Inst%VolcElv(N),z1)
         978           0 :           zTop = MAX(Inst%VolcCld(N),z1)
         979             :
         980             :           ! If volcano is eruptive, zBot /= zTop. In this case, evenly
         981             :           ! distribute emissions in top 1/3 of the plume
         982           0 :           IF ( zBot /= zTop ) THEN
         983           0 :              zBot  = zTop - ( ( zTop - zBot ) / 3.0_sp )
         984           0 :              Erupt = .TRUE.
         985             :           ELSE
         986             :              Erupt = .FALSE.
         987             :           ENDIF
         988             :
         989             :           ! Volcano plume height
         990           0 :           PlumeHgt = zTop - zBot
         991             :
         992             :           ! Distribute emissions into emission arrays. The volcano plume
         993             :           ! ranges from zBot to zTop.
         994           0 :           DO L = 1, HcoState%NZ
         995             :
         996             :              ! Get top height of this box
         997           0 :              z2 = z1 + HcoState%Grid%BXHEIGHT_M%Val(I,J,L)
         998             :
         999             :              ! Skip if the plume bottom is above this grid box top
        1000           0 :              IF ( zBot >= z2 ) THEN
        1001             :                 z1 = z2
        1002             :                 CYCLE
        1003             :              ENDIF
        1004             :
        1005             :              ! If the plume top is below this grid box bottom, we can exit
        1006             :              ! since there will be no more emissions to distribute.
        1007           0 :              IF ( zTop < z1 ) EXIT
        1008             :
        1009             :              ! If we make it to here, the volcano plume is at least partly
        1010             :              ! within this level. Determine the fraction of the plume that
        1011             :              ! is within heights z1 to z2.
        1012             :
        1013             :              ! Get the bottom and top height of the plume within this layer.
        1014           0 :              tmp1 = MAX(z1,zBot)  ! this layer's plume bottom
        1015           0 :              tmp2 = MIN(z2,zTop)  ! this layer's plume top
        1016             :
        1017             :              ! Special case that zTop is heigher than the highest level: make
        1018             :              ! sure that all emissions are going to be used.
        1019           0 :              IF ( ( L == HcoState%NZ ) .AND. ( zTop > z2 ) ) THEN
        1020           0 :                 tmp2 = zTop
        1021             :              ENDIF
        1022             :
        1023             :              ! Fraction of total plume that is within this layer
        1024           0 :              IF ( PlumeHgt == 0.0_sp ) THEN
        1025             :                 Frac = 1.0_sp
        1026             :              ELSE
        1027           0 :                 Frac = (tmp2-tmp1) / PlumeHgt
        1028             :              ENDIF
        1029             :
        1030             :              ! Distribute emissions
        1031           0 :              IF ( Erupt ) THEN
        1032           0 :                 SO2e(I,J,L) = SO2e(I,J,L) + ( Frac * nSo2 )
        1033             :                 volcE       = volcE &
        1034           0 :                             + ( Frac * nSo2 * HcoState%Grid%AREA_M2%Val(I,J) )
        1035             :              ELSE
        1036           0 :                 SO2d(I,J,L) = SO2d(I,J,L) + ( Frac * nSo2 )
        1037             :                 volcD       = volcD &
        1038           0 :                             + ( Frac * nSo2 * HcoState%Grid%AREA_M2%Val(I,J) )
        1039             :              ENDIF
        1040             :
        1041             :              ! The top height is the new bottom
        1042           0 :              z1 = z2
        1043             :           ENDDO
        1044             :
        1045             :           ! testing
        1046             :           !IF ( HCO_IsVerb( HcoState%Config%Err ) ) THEN
        1047             :           !   WRITE(MSG,*) 'Total eruptive  emissions of volcano ', N, ' [kgS/s]: ', volcE
        1048             :           !   CALL HCO_MSG(HcoState%Config%Err,MSG)
        1049             :           !   WRITE(MSG,*) 'Total degassing emissions of volcano ', N, ' [kgS/s]: ', volcD
        1050             :           !   CALL HCO_MSG(HcoState%Config%Err,MSG)
        1051             :           !ENDIF
        1052             :
        1053             :           ! total
        1054           0 :           totE = totE + volcE
        1055           0 :           totD = totD + volcD
        1056             :
        1057             :        ENDDO
        1058             :     ENDIF
        1059             :
        1060             :     ! verbose
        1061           0 :     IF ( HCO_IsVerb( HcoState%Config%Err ) ) THEN
        1062           0 :        WRITE(MSG,*) 'Total eruptive  emissions [kgS/s]: ', totE
        1063           0 :        CALL HCO_MSG(HcoState%Config%Err,MSG)
        1064           0 :        WRITE(MSG,*) 'Total degassing emissions [kgS/s]: ', totD
        1065           0 :        CALL HCO_MSG(HcoState%Config%Err,MSG)
        1066             :     ENDIF
        1067             :
        1068             :     ! Return w/ success
        1069           0 :     RC = HCO_SUCCESS
        1070             :
        1071             :   END SUBROUTINE EmitVolc
        1072             : !EOC
        1073             : !------------------------------------------------------------------------------
        1074             : !                   Harmonized Emissions Component (HEMCO)                    !
        1075             : !------------------------------------------------------------------------------
        1076             : !BOP
        1077             : !
        1078             : ! !IROUTINE: InstGet
        1079             : !
        1080             : ! !DESCRIPTION: Subroutine InstGet returns a poiner to the desired instance.
        1081             : !\\
        1082             : !\\
        1083             : ! !INTERFACE:
        1084             : !
        1085           0 :   SUBROUTINE InstGet( Instance, Inst, RC, PrevInst )
        1086             : !
        1087             : ! !INPUT PARAMETERS:
        1088             : !
        1089             :     INTEGER                             :: Instance
        1090             :     TYPE(MyInst),     POINTER           :: Inst
        1091             :     INTEGER                             :: RC
        1092             :     TYPE(MyInst),     POINTER, OPTIONAL :: PrevInst
        1093             : !
        1094             : ! !REVISION HISTORY:
        1095             : !  18 Feb 2016 - C. Keller   - Initial version
        1096             : !  See https://github.com/geoschem/hemco for complete history
        1097             : !EOP
        1098             : !------------------------------------------------------------------------------
        1099             : !BOC
        1100             :     TYPE(MyInst),     POINTER    :: PrvInst
        1101             :
        1102             :     !=================================================================
        1103             :     ! InstGet begins here!
        1104             :     !=================================================================
        1105             :
        1106             :     ! Get instance. Also archive previous instance.
        1107           0 :     PrvInst => NULL()
        1108           0 :     Inst    => AllInst
        1109           0 :     DO WHILE ( ASSOCIATED(Inst) )
        1110           0 :        IF ( Inst%Instance == Instance ) EXIT
        1111           0 :        PrvInst => Inst
        1112           0 :        Inst    => Inst%NextInst
        1113             :     END DO
        1114           0 :     IF ( .NOT. ASSOCIATED( Inst ) ) THEN
        1115           0 :        RC = HCO_FAIL
        1116           0 :        RETURN
        1117             :     ENDIF
        1118             :
        1119             :     ! Pass output arguments
        1120           0 :     IF ( PRESENT(PrevInst) ) PrevInst => PrvInst
        1121             :
        1122             :     ! Cleanup & Return
        1123           0 :     PrvInst => NULL()
        1124           0 :     RC = HCO_SUCCESS
        1125             :
        1126             :   END SUBROUTINE InstGet
        1127             : !EOC
        1128             : !------------------------------------------------------------------------------
        1129             : !                   Harmonized Emissions Component (HEMCO)                    !
        1130             : !------------------------------------------------------------------------------
        1131             : !BOP
        1132             : !
        1133             : ! !IROUTINE: InstCreate
        1134             : !
        1135             : ! !DESCRIPTION: Subroutine InstCreate creates a new instance.
        1136             : !\\
        1137             : !\\
        1138             : ! !INTERFACE:
        1139             : !
        1140           0 :   SUBROUTINE InstCreate( ExtNr, Instance, Inst, RC )
        1141             : !
        1142             : ! !INPUT PARAMETERS:
        1143             : !
        1144             :     INTEGER,       INTENT(IN)       :: ExtNr
        1145             : !
        1146             : ! !OUTPUT PARAMETERS:
        1147             : !
        1148             :     INTEGER,       INTENT(  OUT)    :: Instance
        1149             :     TYPE(MyInst),  POINTER          :: Inst
        1150             : !
        1151             : ! !INPUT/OUTPUT PARAMETERS:
        1152             : !
        1153             :     INTEGER,       INTENT(INOUT)    :: RC
        1154             : !
        1155             : ! !REVISION HISTORY:
        1156             : !  18 Feb 2016 - C. Keller   - Initial version
        1157             : !  See https://github.com/geoschem/hemco for complete history
        1158             : !EOP
        1159             : !------------------------------------------------------------------------------
        1160             : !BOC
        1161             :     TYPE(MyInst), POINTER          :: TmpInst
        1162             :     INTEGER                        :: nnInst
        1163             :
        1164             :     !=================================================================
        1165             :     ! InstCreate begins here!
        1166             :     !=================================================================
        1167             :
        1168             :     ! ----------------------------------------------------------------
        1169             :     ! Generic instance initialization
        1170             :     ! ----------------------------------------------------------------
        1171             :
        1172             :     ! Initialize
        1173           0 :     Inst => NULL()
        1174             :
        1175             :     ! Get number of already existing instances
        1176           0 :     TmpInst => AllInst
        1177           0 :     nnInst = 0
        1178           0 :     DO WHILE ( ASSOCIATED(TmpInst) )
        1179           0 :        nnInst  =  nnInst + 1
        1180           0 :        TmpInst => TmpInst%NextInst
        1181             :     END DO
        1182             :
        1183             :     ! Create new instance
        1184           0 :     ALLOCATE(Inst)
        1185           0 :     Inst%Instance  = nnInst + 1
        1186           0 :     Inst%ExtNr     = ExtNr
        1187           0 :     Inst%YmdOnFile = -1
        1188             :
        1189             :     ! Attach to instance list
        1190           0 :     Inst%NextInst => AllInst
        1191           0 :     AllInst       => Inst
        1192             :
        1193             :     ! Update output instance
        1194           0 :     Instance = Inst%Instance
        1195             :
        1196             :     ! ----------------------------------------------------------------
        1197             :     ! Type specific initialization statements follow below
        1198             :     ! ----------------------------------------------------------------
        1199             :
        1200             :     ! Return w/ success
        1201           0 :     RC = HCO_SUCCESS
        1202             :
        1203           0 :   END SUBROUTINE InstCreate
        1204             : !EOC
        1205             : !------------------------------------------------------------------------------
        1206             : !                   Harmonized Emissions Component (HEMCO)                    !
        1207             : !------------------------------------------------------------------------------
        1208             : !BOP
        1209             : !
        1210             : ! !IROUTINE: InstRemove
        1211             : !
        1212             : ! !DESCRIPTION: Subroutine InstRemove creates a new instance.
        1213             : !\\
        1214             : !\\
        1215             : ! !INTERFACE:
        1216             : !
        1217           0 :   SUBROUTINE InstRemove( Instance )
        1218             : !
        1219             : ! !INPUT PARAMETERS:
        1220             : !
        1221             :     INTEGER :: Instance
        1222             : !
        1223             : ! !REVISION HISTORY:
        1224             : !  18 Feb 2016 - C. Keller   - Initial version
        1225             : !  See https://github.com/geoschem/hemco for complete history
        1226             : !EOP
        1227             : !------------------------------------------------------------------------------
        1228             : !BOC
        1229             :     INTEGER                     :: RC
        1230             :     TYPE(MyInst), POINTER       :: PrevInst
        1231             :     TYPE(MyInst), POINTER       :: Inst
        1232             :
        1233             :     !=================================================================
        1234             :     ! InstRemove begins here!
        1235             :     !=================================================================
        1236             :
        1237             :     ! Init
        1238           0 :     PrevInst => NULL()
        1239           0 :     Inst     => NULL()
        1240             :
        1241             :     ! Get instance. Also archive previous instance.
        1242           0 :     CALL InstGet ( Instance, Inst, RC, PrevInst=PrevInst )
        1243             :
        1244             :     ! Instance-specific deallocation
        1245           0 :     IF ( ASSOCIATED(Inst) ) THEN
        1246             :
        1247             :        !---------------------------------------------------------------------
        1248             :        ! Deallocate fields of Inst before popping off from the list
        1249             :        ! in order to avoid memory leaks (Bob Yantosca (17 Aug 2022)
        1250             :        !---------------------------------------------------------------------
        1251           0 :        IF ( ALLOCATED( Inst%VolcSlf ) ) THEN
        1252           0 :           DEALLOCATE( Inst%VolcSlf )
        1253             :        ENDIF
        1254             :
        1255           0 :        IF ( ALLOCATED( Inst%VolcElv ) ) THEN
        1256           0 :           DEALLOCATE( Inst%VolcElv )
        1257             :        ENDIF
        1258             :
        1259           0 :        IF ( ALLOCATED( Inst%VolcCld ) ) THEN
        1260           0 :           DEALLOCATE( Inst%VolcCld )
        1261             :        ENDIF
        1262             :
        1263           0 :        IF ( ALLOCATED( Inst%VolcIdx ) ) THEN
        1264           0 :           DEALLOCATE( Inst%VolcIdx )
        1265             :        ENDIF
        1266             :
        1267           0 :        IF ( ALLOCATED( Inst%VolcJdx ) ) THEN
        1268           0 :           DEALLOCATE( Inst%VolcJdx )
        1269             :        ENDIF
        1270             :
        1271           0 :        IF ( ALLOCATED( Inst%SpcIDs ) ) THEN
        1272           0 :           DEALLOCATE( Inst%SpcIDs )
        1273             :        ENDIF
        1274             :
        1275           0 :        IF ( ALLOCATED( Inst%SpcScl ) ) THEN
        1276           0 :           DEALLOCATE( Inst%SpcScl )
        1277             :        ENDIF
        1278             :
        1279           0 :        IF ( ALLOCATED( Inst%SpcScalFldNme ) ) THEN
        1280           0 :           DEALLOCATE( Inst%SpcScalFldNme )
        1281             :        ENDIF
        1282             :
        1283             :        !---------------------------------------------------------------------
        1284             :        ! Pop off instance from list
        1285             :        !---------------------------------------------------------------------
        1286           0 :        IF ( ASSOCIATED(PrevInst) ) THEN
        1287           0 :           PrevInst%NextInst => Inst%NextInst
        1288             :        ELSE
        1289           0 :           AllInst => Inst%NextInst
        1290             :        ENDIF
        1291           0 :        DEALLOCATE(Inst)
        1292             :     ENDIF
        1293             :
        1294             :     ! Free pointers before exiting
        1295           0 :     PrevInst => NULL()
        1296           0 :     Inst     => NULL()
        1297             :
        1298           0 :    END SUBROUTINE InstRemove
        1299             : !EOC
        1300           0 : END MODULE HCOX_Volcano_Mod

Generated by: LCOV version 1.14