# Creating County Mask for DESID using shp2cmaq

**URL:** <https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499>\
**Category:** DESID\
**Created:** [January 7, 2025, 7:42pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499 "2025-01-07T19:42:20Z")\
**Posts on this page:** 13\
**Page:** 1

<div class="post-metadata">

**Author:** ![DEEP](https://avatars.discourse-cdn.com/v4/letter/d/ce7236/32.png) [@DEEP](https://forum.cmascenter.org/u/DEEP)\
**Post date:** [January 7, 2025, 7:42pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/1 "2025-01-07T19:42:20Z")

</div>

I am trying to use shp2cmaq to create a county mask file to use with DESID. I would then use DESID to reduce emissions in a group of counties that would be a region family. I am trying to create the mask for the CMAQ version 5.5 benchmark which uses the 2018\_12NE3 grid. The county shape file I am using comes from a US Census tigerline file:

[https://www.census.gov/cgi-bin/geo/shapefiles/index.php?year=2020&layergroup=Counties+%28and+equivalent%29](https://www.census.gov/cgi-bin/geo/shapefiles/index.php?year=2020&layergroup=Counties+%28and+equivalent%29)

The tiger file contains the 5 digit state/county fips code in a column called GEOID. I call shp2cmaq using the tigerline zipfile with attrkey=‘GEOID’, gdnam=‘2018\_12NE3’,  
gdpath=‘/…/CMAQv5.4\_2018\_12NE3\_Benchmark\_2Day\_Input/2018\_12NE3/GRIDDESC’,  
outpath=‘/…/CountyMask.nc’,  
outformat=‘NETCDF3\_CLASSIC’, srckey=‘intersection\_area’,overwrite=True,verbose=1. The crs for the file is EPSG:4269.

The CountyMask.nc file that is created does not seem to function properly. Examining it with Panoply shows that the counties are overlapped with the grid in unknown units. Panoply output:

netcdf /home/…/CountyMask.nc {  
dimensions:  
TSTEP = 1;  
VAR = 667;  
DATE-TIME = 2;  
LAY = 1;  
ROW = 105;  
COL = 100;  
variables:  
int TFLAG(TSTEP=1, VAR=667, DATE-TIME=2);  
:units = “\<YYYYJJJ,HHMMSS\>”;  
:long\_name = "TFLAG ";  
:var\_desc = "TFLAG ";

```
float GEOID_09001(TSTEP=1, LAY=1, ROW=105, COL=100);
  :_FillValue = NaNf; // float
  :long_name = "GEOID_09001 ";
  :var_desc = "Fractional overlap of GEOID_09001 with 2018_12NE3 ";
  :unit = "1 ";
  :units = "unknown";

float GEOID_09003(TSTEP=1, LAY=1, ROW=105, COL=100);
  :_FillValue = NaNf; // float
  :long_name = "GEOID_09003 ";
  :var_desc = "Fractional overlap of GEOID_09003 with 2018_12NE3 ";
  :unit = "1 ";
  :units = "unknown";

```

Whereas Panoply indicates that the gridmask for states that comes with the Benchmark has the units of “fraction” as below:

netcdf /data1/CMAQ/DATA/CMAQv5.4\_2018\_12NE3\_Benchmark\_2Day\_Input/2018\_12NE3/GRIDMASK\_STATES\_12NE3.nc {  
dimensions:  
TSTEP = 1;  
DATE-TIME = 2;  
LAY = 1;  
VAR = 49;  
ROW = 105;  
COL = 100;  
variables:  
int TFLAG(TSTEP=1, VAR=49, DATE-TIME=2);  
:units = “\<YYYYDDD,HHMMSS\>”;  
:long\_name = "TFLAG ";  
:var\_desc = "Timestep-valid flags: (1) YYYYDDD or (2) HHMMSS ";

```
float AL(TSTEP=1, LAY=1, ROW=105, COL=100);
  :long_name = "AL ";
  :units = "fraction";
  :var_desc = "AL fractional area per grid cell ";

float AZ(TSTEP=1, LAY=1, ROW=105, COL=100);
  :long_name = "AZ ";
  :units = "fraction";
  :var_desc = "AZ fractional area per grid cell ";

```

How am I misusing shp2cmaq?  
Thank you.

---

<div class="post-metadata">

**Author:** ![barronh](https://avatars.discourse-cdn.com/v4/letter/b/85e7bf/32.png) [@barronh](https://forum.cmascenter.org/u/barronh)\
**Post date:** [January 16, 2025, 4:01pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/2 "2025-01-16T16:01:37Z")

</div>

Thanks for the question and posting it. You’re right that unit/units should be fixed, but I don’t think that’s driving the problem.

For the most part, I don’t think that var\_desc or units are used by DESID. Thought @Ben_Murphy might correct me. So, I think that it is something else. If you want to fix that, edit the shp2cmaq.py file to change unit to units on line 233ish. That will overwrite unknown with “1”, which is an appropriate unit. If you really want it to be fraction, you can change it from 1 to fraction. Again, I don’t think this will help – but easy to test.

I’m glad you’re using the NETCDF3\_CLASSIC, which rules out the problem where the netcdf library CMAQ uses is not compatible. A few more questions to help figure this out:

- What does it mean that “CountyMask.nc … does not seem to function properly”? Do you get an error? Or does CMAQ not apply the changes you expect?
- It might be useful to see your DESID config and the CMAQ log.

For the default file, I added total and dominant variables with the code below:

```auto
keys = sorted([k for k in f if k != 'TFLAG' and not k.startswith('STATE_')])
f['STATE_TOT'] = eval(' + '.join([k for k in keys if k != 'TFLAG']), None, f.data_vars)
f['STATE_DOM'] = f['STATE_TOT'].dims, np.argmax(np.stack([f[k] for k in keys]), axis=0)

```

And plotted:

 ![image](https://canada1.discourse-cdn.com/flex027/uploads/cmas/original/2X/5/5cd90d1c79fb56e45ac29c085c90797a9577a295.png)

And compared to what you shared to me via email:

 ![county](https://canada1.discourse-cdn.com/flex027/uploads/cmas/original/2X/3/37dda3ae749df8d0dd08f9e75b787636e3919d21.jpeg)

The domains generally match with expected coverages. So, it looks like a lot is going right…

---

<div class="post-metadata">

**Author:** ![DEEP](https://avatars.discourse-cdn.com/v4/letter/d/ce7236/32.png) [@DEEP](https://forum.cmascenter.org/u/DEEP)\
**Post date:** [January 16, 2025, 9:17pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/3 "2025-01-16T21:17:22Z")

</div>

CMAQ does not apply the changes I expect. In Verdi the NO and NO2 that I had intended to zero out are identical to the benchmark reference case and are clearly not zeroed out.

This is the runscript file:  
[run\_Bench\_CountyMask.csh](https://forum.cmascenter.org/uploads/short-url/3w91Xp9JmSrOeunRafCoimdyfnQ.csh) (39.5 KB

Note that for the following files I added .txt to end of the filename so they would upload to this forum.

This is a log file:  
CTM\_LOG\_000.v55\_gcc\_Bench\_2018\_12NE3\_cb6r5\_ae7\_aq\_m3dry\_20180701  
[CTM\_LOG\_000.v55\_gcc\_Bench\_2018\_12NE3\_cb6r5\_ae7\_aq\_m3dry\_20180701.txt](https://forum.cmascenter.org/uploads/short-url/s3WOvzDZ8HZOXjjq2aHZtIKlN0K.txt) (709.1 KB)

These are the nml files that describe the region family and zero the NOx emissions:  
[CMAQ\_Control\_DESID.txt](https://forum.cmascenter.org/uploads/short-url/hpgGxpAzmXJvnmrEW3WVA0XY6R.txt) (15.6 KB)

and  
[CMAQ\_Control\_DESID\_cb6r5\_ae7\_aq.txt](https://forum.cmascenter.org/uploads/short-url/tYKTPQtXow2LmxFiPJ6uciNpYBL.txt) (13.2 KB)

This is the run log:  
[cctmCOUNTYlog.txt](https://forum.cmascenter.org/uploads/short-url/7W90c3DxJj02ZrJgjse6lviiJd6.txt) (493.1 KB)

Thank you.

---

<div class="post-metadata">

**Author:** ![lizadams](https://avatars.discourse-cdn.com/v4/letter/l/49beb7/32.png) [@lizadams](https://forum.cmascenter.org/u/lizadams)\
**Post date:** [January 16, 2025, 11:01pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/4 "2025-01-16T23:01:56Z")

</div>

This is the file you are using in your run script:

#\> Spatial Masks For Emissions Scaling  
setenv CMAQ\_MASKS $INPDIR/GRIDMASK\_COUNTIES\_12NE3.nc

I assume this is the file that you generated using following these instructions: [CMAQ/DOCS/Users\_Guide/Tutorials/CMAQ\_UG\_tutorial\_emissions.md at main · USEPA/CMAQ · GitHub](https://github.com/USEPA/CMAQ/blob/main/DOCS/Users_Guide/Tutorials/CMAQ_UG_tutorial_emissions.md#11-apply-scaling-with-spatial-dependence)

Using the following tool:

> **[CMAQ/PYTOOLS/shp2cmaq at main · USEPA/CMAQ](https://github.com/USEPA/CMAQ/tree/main/PYTOOLS/shp2cmaq)**
>
> Code for U.S. EPA’s Community Multiscale Air Quality Model (CMAQ) which helps in conducting air quality model simulations - USEPA/CMAQ

Your CTM\_LOG file contains the following WARNINGs

```auto
     |> Reading Emission Control Namelist:
     +======================================

     Performing Basic Error Checks for Emission Scaling Rules

     "CMAQ_MASKS" opened as OLD:READ-ONLY
     File name "/data1/CMAQ/DATA/CMAQv5.4_2018_12NE3_Benchmark_2Day_Input/2018_12NE3/GRIDMASK_COUNTIES_12NE3.nc"
     File type GRDDED3
     Execution ID "NA"
     Grid name "12NE3"
     Dimensions: 105 rows, 100 cols, 1 lays, 667 vbles
     NetCDF ID: 524288 opened as READONLY
     Time-independent data.

 >>--->> WARNING in subroutine RDTFLAG
     Time step not available in file CMAQ_MASKS for variable GEOID_09001

     >>--->> WARNING in subroutine XTRACT3
     Time step not available for file: CMAQ_MASKS

     >>--->> WARNING in subroutine DESID_INIT_REGIO on PE 000
     Could not read GEOID_09001 from file CMAQ_MASKS

     >>--->> WARNING in subroutine RDTFLAG
     Time step not available in file CMAQ_MASKS for variable GEOID_09007

... this warning was repeated for all of the defined regions GEOID_*

     >>--->> WARNING in subroutine DESID_INIT_REGIO on PE 000
     Could not read GEOID_34041 from file CMAQ_MASKS

     Closing file CMAQ_MASKS

```

I am not sure why the GEOID variables can’t be read by CMAQ, but if the regions can’t be read from the CMAQ\_MASKS file, then the DESID won’t be able to apply the emission reduction.

One thing you can try is to use m3xtract to recreate the file GRIDMASK\_COUNTIES\_12NE3.nc file.

setenv INFILE GRIDMASK\_COUNTIES\_12NE3.nc  
setenv OUTFILE GRIDMASK\_COUNTIES\_12NE3\_redo.nc

m3xtract  
(extract all layers, all variables)

---

<div class="post-metadata">

**Author:** ![hogrefe.christian](https://avatars.discourse-cdn.com/v4/letter/h/74df32/32.png) [@hogrefe.christian](https://forum.cmascenter.org/u/hogrefe.christian)\
**Post date:** [January 17, 2025, 2:13am UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/5 "2025-01-17T02:13:56Z")

</div>

To better understand why CMAQ fails to read the mask file, could you please upload a text file showing the output of the command `ncdump -v TFLAG GRIDMASK_COUNTIES_12NE3.nc`?

---

<div class="post-metadata">

**Author:** ![DEEP](https://avatars.discourse-cdn.com/v4/letter/d/ce7236/32.png) [@DEEP](https://forum.cmascenter.org/u/DEEP)\
**Post date:** [January 17, 2025, 2:09pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/6 "2025-01-17T14:09:12Z")

</div>

Yes, I have read the DESID and shp2cmaq documentation on the github page. I noticed the “could not read” warnings and was concerned that my problem was the variable name as shp2cmaq indicates that variable length and certain characters could cause problems. I tried several things including removing the “GEOID\_” prefix, but that didn’t help and I was concerned that the string name might then be interpreted as numeric and the leading zeros of the 5 digit FIPS code might be lost. I eventually concluded that the variable name was not the problem and began focusing on the units “unknown”. I also conducted tests to create a region family with the benchmark state gridmask to assure I was understanding that properly. I then tried creating my own state gridmask with shp2cmaq and the file listed in the shp2cmaq instructions (cb\_2022\_us\_state\_500k.zip). The results of that were similar to the county results.

I should explain that I create CountyMask.nc in my home directory and copy it to GRIDMASK\_COUNTIES\_12NE3.nc in the benchmark input directory.

I will look into using m3xtract.

Thank you.

---

<div class="post-metadata">

**Author:** ![DEEP](https://avatars.discourse-cdn.com/v4/letter/d/ce7236/32.png) [@DEEP](https://forum.cmascenter.org/u/DEEP)\
**Post date:** [January 17, 2025, 2:12pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/7 "2025-01-17T14:12:24Z")

</div>

Attached is the out from ncdump -v TFLAG GRIDMASK\_COUNTIES\_12NE3.nc

[ncdump.txt](https://forum.cmascenter.org/uploads/short-url/vK1oVylD3SzB7WQ4RxCfmf3WMy4.txt) (241.8 KB)

Thank you.

---

<div class="post-metadata">

**Author:** ![barronh](https://avatars.discourse-cdn.com/v4/letter/b/85e7bf/32.png) [@barronh](https://forum.cmascenter.org/u/barronh)\
**Post date:** [January 17, 2025, 4:06pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/8 "2025-01-17T16:06:25Z")

</div>

I think the problem is probably related to the time-independence notation. You’ll need to update the gpath to point to your file. I believe the file is not being identified as time-independent. The code below should work in your python environment and I think it will fix the file for you.

```python
import netCDF4
import numpy as np
gpath = 'GRIDMASK_COUNTIES_12NE3.nc'
f = netCDF4.Dataset(gpath, mode='a')
f.SDATE = np.int32(-635)
f.variables['TFLAG'][:, :, 0] = np.int32(0)
f.close()

```

After you confirm that this works, I’ll dig a bit deeper because the raw output seems to be working for others.

---

<div class="post-metadata">

**Author:** ![barronh](https://avatars.discourse-cdn.com/v4/letter/b/85e7bf/32.png) [@barronh](https://forum.cmascenter.org/u/barronh)\
**Post date:** [January 17, 2025, 4:43pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/9 "2025-01-17T16:43:11Z")

</div>

Can you confirm that you’re using the version of shp2cmaq that has a shp2cmaq.py and not just the Notebook (.ipynb)? And, please confirm the version of cmaqsatproc that you’re using.

```python
import cmaqsatproc as csp
print('cmaqsatproc', csp. __version__ )
import shp2cmaq
print('shp2cmaq', shp2cmaq. __version__ )

```

For example, when I run this code, I get:

```auto
cmaqsatproc 0.4.1
shp2cmaq 1.0

```

---

<div class="post-metadata">

**Author:** ![DEEP](https://avatars.discourse-cdn.com/v4/letter/d/ce7236/32.png) [@DEEP](https://forum.cmascenter.org/u/DEEP)\
**Post date:** [January 17, 2025, 5:44pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/10 "2025-01-17T17:44:08Z")

</div>

Confirming that versions are correct:

 ![](https://canada1.discourse-cdn.com/flex027/uploads/cmas/original/2X/8/89b7abb1b7011426dd418489418f0a895529cdde.png)

---

<div class="post-metadata">

**Author:** ![DEEP](https://avatars.discourse-cdn.com/v4/letter/d/ce7236/32.png) [@DEEP](https://forum.cmascenter.org/u/DEEP)\
**Post date:** [January 23, 2025, 8:35pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/11 "2025-01-23T20:35:27Z")

</div>

This solved it. Using m3xtract to essentially copy the file created by shp2cmaq allowed the county shapes (GEIOD\_xxxxx) to be read from the CMAQ\_MASKS file.

Thank you all for your help and patience.

---

<div class="post-metadata">

**Author:** ![DEEP](https://avatars.discourse-cdn.com/v4/letter/d/ce7236/32.png) [@DEEP](https://forum.cmascenter.org/u/DEEP)\
**Post date:** [January 23, 2025, 8:42pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/12 "2025-01-23T20:42:11Z")

</div>

I did also use this in the solution. I ran this code (I used mode =‘r+’) after creating the mask with shp2cmaq, and then ran that file through m3xtract.  
Thanks.

---

<div class="post-metadata">

**Author:** ![DEEP](https://avatars.discourse-cdn.com/v4/letter/d/ce7236/32.png) [@DEEP](https://forum.cmascenter.org/u/DEEP)\
**Post date:** [January 31, 2025, 3:20pm UTC](https://forum.cmascenter.org/t/creating-county-mask-for-desid-using-shp2cmaq/5499/13 "2025-01-31T15:20:14Z")

</div>

Having obtained the correct results by applying the timestep and m3xtract fixes, I wanted to see if either alone provided the correction. In doing so I could not replicate my original error. Upon further investigation and recalling that I used the netcdf3 option in shp2cmaq after a warning about HDF5, I found that HDF5 files had been updated just prior to my obtaining the correct results. I therefore suspect the update to HDF5 may have been a contributing factor to obtaining the correct results.
