summation of netcdf grids

How do I...?
Post Reply
User avatar
kbedka1
Posts: 428
Joined: Wed Jan 28, 2009 7:27 pm

summation of netcdf grids

Post by kbedka1 »

I have 5 gridded netcdf files with GLM lightning that have no timestamps or time coordinate so I can't aggregate like I normally would with existing Mc-V formulas. Do you have a recommendation for how I can simply sum these 5 files to get a 5-min accumulation of GLM lightning flashes? My scripting ability is poor so any recipes for this process would be helpful. I've attached the 5 files I mention
Attachments
IXTR99_KNES_021514_8788.2018080215.nc
(117.33 KiB) Downloaded 272 times
IXTR99_KNES_021515_8790.2018080215.nc
(117.2 KiB) Downloaded 286 times
IXTR99_KNES_021516_8792.2018080215.nc
(117.36 KiB) Downloaded 289 times
IXTR99_KNES_021517_8794.2018080215.nc
(117.55 KiB) Downloaded 272 times
IXTR99_KNES_021518_8796.2018080215.nc
(117.17 KiB) Downloaded 294 times
User avatar
bobc
Posts: 988
Joined: Mon Nov 15, 2010 5:57 pm

Re: summation of netcdf grids

Post by bobc »

Hi Kris,

The easiest way to do this with things as they are is to run a script through the Jython Shell. You can use this script as an example. All you should have to do is change the dataDir variable to point to the data directory on your machine. Note that all of the *.nc files in the dataDir directory will be loaded and displayed:

Code: Select all

# import glob, which will be used later
import glob

# define the directory where the data exists
dataDir = 'C:/Users/rcarp/Data/GLM/'

# initialize an empty list.  data objects will be added to this list later
dataList = []

# loop through all of the *.nc files in the dataDir directory and load
# the Flash_extent_density field of each of them.  Append each data object
# to the dataList list initialized earlier.
for f in glob.glob('%s*.nc' % dataDir):
    dataParms = dict(
        filename = f,
        field = 'Flash_extent_density'
        )
    # loadGrid function for each file
    data = loadGrid(**dataParms)
    # append each data object to the dataList list
    dataList.append(data)

# display a loop of the dataList
layer = activeDisplay().createLayer('Image Sequence Display', dataList)

Normally, when the filename includes the date/time information you can write a NcML wrapper file where you can use dateFormatMark to extract the date/time information. This allows you to simply go to the Data Sources tab of the Data Explorer, select the NcML file, click Add Source, and go on and display your data without a script. However, this doesn't work with these particular files because of the characters separating the date/time. For example:

IXTR99_KNES_021514_8788.2018080215.nc
IXTR99_KNES_021515_8790.2018080215.nc
IXTR99_KNES_021516_8792.2018080215.nc
IXTR99_KNES_021517_8794.2018080215.nc
IXTR99_KNES_021518_8796.2018080215.nc

The red text above contains the hours/minutes. The blue text above contains the year, month, day, and hour. Because the end of each filename string is the same, the minutes would need to be specified in the dateFormatMark. However, this doesn't work because of the 8788, 8790, 8792, 8794, and 8796 parts of the filename string. If you have the ability to change these filenames, it could be worth appending the minutes to the end of the filename strings. This would make the NcML wrapper file with dateFormatMark possible. For example, your files could be named:

IXTR99_KNES_021514_8788.201808021514.nc
IXTR99_KNES_021515_8790.201808021515.nc
IXTR99_KNES_021516_8792.201808021516.nc
IXTR99_KNES_021517_8794.201808021517.nc
IXTR99_KNES_021518_8796.201808021518.nc

You can find more information about dateFormatMark by searching the forums. For example:
McIDAS-V version 1.7u1 and GOES-R
Create time Series from date/time in filename (ncml wrapper)
problem making time sequence, single time grids

Please let me know if you have any questions about this.

Thanks,
Bob Carp
User avatar
kbedka1
Posts: 428
Joined: Wed Jan 28, 2009 7:27 pm

Re: summation of netcdf grids

Post by kbedka1 »

getting back to this, when the data is summed say over 5 timesteps, and a fill value is encountered in 1 of the steps, is the summed count set to be fill or does the summing ignore the fill value? How would the data be summed using the script above?
User avatar
joleenf
Posts: 1123
Joined: Mon Jan 19, 2009 7:16 pm

Re: summation of netcdf grids

Post by joleenf »

Hi Kris,

The script above is an image sequence, not a sum. If these values are summed, missing values will result in a missing value for the grid cell. Were you looking for an image sequence or a single image representing the sum of the five files. My understanding was that you were looking for the latter.

Joleen
User avatar
bobc
Posts: 988
Joined: Mon Nov 15, 2010 5:57 pm

Re: summation of netcdf grids

Post by bobc »

You should be able to use the sum() function to add up the data values at each timestep. For example, at the end of the script, you could add the following two lines:

Code: Select all

dataSum = sum((dataList[0],dataList[1],dataList[2],dataList[3],dataList[4]))
sumLayer = activeDisplay().createLayer('Image Display', dataSum)

I verified that the values from this sum layer are correct at several pixels when comparing to the individual timesteps at the same pixel.

The _FillValue attribute of the Flash_extent_density variable is set to 0UB in the files attached to this post earlier. It appears that the summed layer contains all of the data values at each timestep... meaning that if one timestep has a fill value at a pixel, but another contains real data, then real data will be included in the summed layer.

Does this make sense?

Bob
User avatar
kbedka1
Posts: 428
Joined: Wed Jan 28, 2009 7:27 pm

Re: summation of netcdf grids

Post by kbedka1 »

This made sense and worked for me. How could I then apply the parallax shift formula to the summed array so that it aligns with radar imagery? Also I'd want a contour option. Can I add this somehow to my jython library so I have more functionality with how it works?
User avatar
bobc
Posts: 988
Joined: Mon Nov 15, 2010 5:57 pm

Re: summation of netcdf grids

Post by bobc »

Hi Kris,

I found a way to get this working through the GUI, not using scripting at all. The first thing you'll need to do is copy/paste this text into an NcML file.

Code: Select all

<netcdf xmlns="http://www.unidata.ucar.edu/namespaces/netcdf/ncml-2.2">
  <aggregation dimName="time" type="joinNew">
    <variableAgg name="Flash_extent_density_window"/>
    <variableAgg name="DQF"/>
    <variableAgg name="Flash_extent_density"/>
    <scan location="" dateFormatMark="IXTR99_KNES_#ddHHmm_wwSS.yyyyMM" suffix=".nc" subdirs="false"/>
  </aggregation>
</netcdf>

You could name this whatever you want, as long as it ends with ".ncml". For example, "IXTR99.ncml". This file uses "dateFormatMark" to extract date/time info from the file names. As noted in a previous post on this thread (April 05, 2019), these filenames don't contain all of the date/time characters right next to each other. The biggest issue is the 4-character string before the period (e.g. "8788"). I had to fake things a bit by replacing these characters with "wwSS" for "week" and "millisecond". These "wwSS" characters don't end up changing the timestamp values. If you put this NcML file in the same directory as your netCDF files, you can go to the General>Files/Directories chooser, select the NcML file, and click Add Source. This gets loaded in as a gridded data source.

When you mention parallax correction, I'm assuming that you're referring to the lat/lon shift functions in this forum post. Since we're aggregating netCDF files, I'm using the nc_latlonshift() function and formula, which are in the forum I linked to. Now that your NcML data source is added to the Field Selector, you can do the following:

  1. From the Field Selector, select "Formulas" on the left side of the window under Data Sources.
  2. Choose the "Grids > Time Steps > Total sum over time steps" formula.
  3. Choose the display type you want to use. I used Color-Shaded Plan View.
  4. Click Create Display.
  5. In the Field Selector window that appears, choose the nc_latlonshift formula and click OK.
  6. In the Select Input window, enter your values for latshift and lonshift. Click OK.
  7. In the Field Selector window that appears, expand the NcML data source and select "GLM Flash Extent Density". Click OK.

To verify that the shifted display was correct, I displayed the original summed data without any shift. This was done by choosing the "Grids > Time Steps > Total sum over time steps" formula, and the "GLM Flash Extent Density" field. The shifted data was displaced by the lat/lon I specified when going through the steps above.

How does this work for you?

Bob
Post Reply