Applying a mask to a scale factor does not work

[msulprizio]

Your name

Melissa Sulprizio

Your affiliation

Harvard

What happened? What did you expect to happen?

The example for scaling or zeroing emissions with a rectangular mask does not work. The current example in that text is:

In Base emissions

0 HTAP_NO_IND /path/to/HTAP_NO_INDUSTRY.generic.01x01.nc emi_no 2008-2010/1-12/1/0 C xy kg/m2/s NO 1/27/25/501 1/2 4

In Scale Factors:

501 SCALE_AUS 0.0 - - - xy unitless 1 1010

In Masks:

# Defines a rectangular region that should cover AUS + surrounding islands
1010 AUS_MASK 105.0/-46.0/160.0/-10.0 – 2000/1/1/0 C xy 1 1 105.0/-46.0/160.0/-10.0

I followed the example and modified the following lines in HEMCO_Config.rc to only apply the OH_pert_factor over the Southern Hemisphere:

## Field
* GLOBAL_OH  $ROOT/OH/v2014-09/v5-07-08/OH_3Dglobal.geos5.72L.4x5.nc OH   1985/1-12/1/0 C xyz kg/m3 * 2 1 1

## Scale factor
2 OH_pert_factor  1.1 - - - xy 1 1020

## Mask
1020 S_HEMIS -180/-90/180/0 - 2000/1/1/0 C xy 1 1 -180/-90/180/0

However, I get the error:

HEMCO: Entering GET_CURRENT_EMISSIONS (hco_calc_mod.F90) ( 1)
                                                                               
 Evaluate field GLOBAL_OH
Applying scale factor OH_pert_factor
 
HEMCO ERROR:  Illegal operator for: OH_pert_factor operation:         1020
 --> LOCATION: Apply_Scale_Factor (src/Core/hco_calc_mod.F90)

HEMCO ERROR: Illegal mathematical operator for scale factor:OH_pert_factor
 --> LOCATION: GET_CURRENT_EMISSIONS (hco_calc_mod.F90)
 
HEMCO ERROR: ERROR 14
 --> LOCATION: HCO_EvalFld_3d (HCO_calc_mod.F90)
===============================================================================
GEOS-Chem ERROR: GLOBAL_OH not found in HEMCO data list!
 -> at CHEMCH4 (in module GeosCore/global_ch4_mod.F90)
===============================================================================

===============================================================================
GEOS-Chem ERROR: Error encountered in "ChemCh4"!
 -> at Do_Chemistry  (in module GeosCore/chemistry_mod.F90)
===============================================================================

===============================================================================
GEOS-CHEM ERROR: Error encountered in "Do_Chemistry"!
STOP at  -> at GEOS-Chem (in GeosCore/main.F90)
===============================================================================
     - CLEANUP: deallocating arrays now...

-------------------------------------------------------------------------------
HEMCO 3.9.1 FINISHED.
      Warnings: *
-------------------------------------------------------------------------------

What are the steps to reproduce the bug?

See above

Please attach any relevant configuration and log files.

No response

What HEMCO version were you using?

3.9.1

What environment were you running HEMCO on?

Local cluster

What compiler and version were you using?

gcc 12.2.0

Will you be addressing this bug yourself?

Yes

In what configuration were you running HEMCO?

GCClassic

As what resolution were you running HEMCO?

2x2.5

What meterology fields did you use?

GEOS-FP

Additional information

In routine APPLY_SCALE_FACTOR in hco_calc_mod.F90 there is a check for the operator (i.e. the last integer defined in the scale factor line in HEMCO_Config.rc). In the example case in the docs it tells users to put the mask number there, however the code appears to only allow for three operators (multiply, divide, square).

HEMCO/src/Core/hco_calc_mod.F90

Lines 1585 to 1674 in e68ebae

	!------------------------------------------------------------------------------
	! GEOS-Chem Global Chemical Transport Model !
	!------------------------------------------------------------------------------
	!BOP
	!
	! !IROUTINE: Apply_Scale_Factor
	!
	! !DESCRIPTION: Applies scale factors to the input
	!\\
	!\\
	! !INTERFACE:
	!
	SUBROUTINE Apply_Scale_Factor( I, J, L, ScalDct, scalFac, dataVal, RC )
	!
	! !USES:
	!
	USE HCO_Error_Mod
	!
	! !INPUT PARAMETERS:
	!
	INTEGER, INTENT(IN) :: I ! Longitude (or X-dim) index
	INTEGER, INTENT(IN) :: J ! Latitude (or Y-dim) index
	INTEGER, INTENT(IN) :: L ! Vertical level index
	TYPE(DataCont), POINTER :: ScalDct ! Scale Factor data container
	REAL(sp), INTENT(IN) :: scalFac ! Scale factor
	!
	! !INPUT/OUTPUT PARAMETERS:
	!
	REAL(hp), INTENT(INOUT) :: dataVal ! Data to be scaled
	!
	! !OUTPUT PARAMETERS:
	!
	INTEGER, INTENT(OUT) :: RC ! Success or failure
	!
	! !REMARKS:
	! This code was abstracted out of Get_Current_Emisssions for clarity. We
	! have also refactored the code to remove ELSE blocks for better efficiency.
	!EOP
	!------------------------------------------------------------------------------
	!BOC
	!
	! !LOCAL VARIABLES:
	!
	! Strings
	CHARACTER(LEN=255) :: errMsg, thisLoc
	!========================================================================
	! Apply_Scale_Factor begins here!
	!========================================================================
	! Initialize
	RC = HCO_SUCCESS
	errMsg = ''
	thisLoc = 'Apply_Scale_Factor (src/Core/hco_calc_mod.F90)'
	!------------------------------------------------------------------------
	! Operation = 1: multiply
	!------------------------------------------------------------------------
	IF ( ScalDct%Oper == 1 ) THEN
	dataVal = dataVal * scalFac
	RETURN
	ENDIF
	!------------------------------------------------------------------------
	! Operation = -1: divide
	!------------------------------------------------------------------------
	IF ( ScalDct%Oper == -1 ) THEN
	IF ( scalFac /= 0.0_sp ) THEN
	dataVal = dataVal / scalFac
	ENDIF
	RETURN
	ENDIF
	!------------------------------------------------------------------------
	! Operation = 2: square
	!------------------------------------------------------------------------
	IF ( ScalDct%Oper == 2 ) THEN
	dataVal = dataVal * scalFac * scalFac
	RETURN
	ENDIF
	!------------------------------------------------------------------------
	! Return w/ error otherwise (Oper 3 is only allowed for masks!)
	!------------------------------------------------------------------------
	WRITE( errMsg, * ) 'Illegal operator for: ', TRIM( ScalDct%cName ), &
	' operation: ', ScalDct%Oper
	CALL HCO_Error( ErrMsg, RC, thisLoc )
	RC = 2
	END SUBROUTINE Apply_Scale_Factor

This leads me to believe that masking scale factors is not possible after all.

[yantosca]

Tagging @ktravis213 @jimmielin

[msulprizio]

I was able to get around this issue by creating a hemispheric mask file with the northern hemisphere defined as ID=1 and the southern hemisphere defined as ID=2. I then applied the scaling following the method described in Scale (or zero) emissions with a shapefile country mask.

This is a sufficient workaround for the IMI but perhaps we can consider expanding the code in hco_calc_mod.F90 to allow for masking of scale factors. In the meantime, we should probably remove this section from the HEMCO examples documentation since it doesn't work:

• https://hemco.readthedocs.io/en/stable/hco-ref-guide/more-examples.html#scale-or-zero-emissions-with-a-rectangular-mask

[yantosca]

Thanks @msulprizio, I will remove that text.

[yantosca]

@msulprizio: The text has been removed from RTD in commit a657308 (currently in the docs/dev branch).

[msulprizio]

Thanks @yantosca!

[github-actions]

This issue has been automatically marked as stale because it has not had recent activity. If there are no updates within 7 days it will be closed. You can add the "never stale" tag to prevent the issue from closing this issue.