[ncl-talk] Problems converting hdf to netcdf using NCL
Dennis Shea
shea at ucar.edu
Wed May 27 13:28:29 MDT 2015
The most important rule in data processing is 'look at your data'.
[1] Swath data often has missing values. Just look at a sample source .hdf
file
%> ncl_filedump -v Optical_Depth_Land_And_Ocean
MOD04_L2.A2014305.0000.006.2015077213911.pscs_000500937810.hdf | less
[snip]
Variable No. 0
Optical_Depth_Land_And_Ocean:
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 3 24 -9999 65 42 -9999 280
257
249 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 105 -9999 -9999 -9999 -9999 127 123
157
167 -9999 117 112 125 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 7
-9999
-9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 53 53
56
[snip
the -9999 are missing values (_FillValue).
If you want the whole globe filled with data ... you will have to fill it
in.
[2] I don't use ncview but I think is is looking at one time step.
The NCL plot script should plot all the available data at once.
The files you sent are not labeled. I have no idea what each refers to.
On Wed, May 27, 2015 at 12:39 PM, BLIUJUS, STEVEN D CTR USAF AFMC
AFLCMC/HBAW-OL <steven.bliujus.3.ctr at us.af.mil> wrote:
> Found the error Dennis, there was a discrepancy in the dimension size. One
> of the files has a 205 dimension, so I reformatted the do loop section to
> look like:
>
> mAcross = 135
> nAlong = 203
> nAlong1 = nAlong + 2
>
> time = new ( nfili , "double", "No_FillValue" )
> data = new ( (/nfili,nAlong1,mAcross/), "float" , -9999.0)
> lat2d = new ( (/nfili,nAlong1,mAcross/), "float" , -999.0)
> lon2d = new ( (/nfili,nAlong1,mAcross/), "float" , -999.0)
> print("nfili = ")
> print(nfili)
> do nf=0,nfili-1 ; loop over files
> a = addfile(diri+fili(nf),"r")
> time(nf) = (/ parse_MODIS( fili(nf) ) /)
> print("nf="+nf+" time(nf)="+time(nf))
>
> x := short2flt(a->$var$)
> dimx = dimsizes(x)
> print("dimsizex="+dimx)
>
> if (dimx(0).eq.203 .and. dimx(1).eq.135) then
> data(nf,0:202,:) = (/ x /)
> lat2d(nf,0:202,:) = (/ a->Latitude /)
> lon2d(nf,0:202,:) = (/ a->Longitude /)
> else if(dimx(0).eq.204 .and. dimx(1).eq.135) then
> data(nf,0:203:,:) = (/ x /)
> lat2d(nf,0:203,:) = (/ a->Latitude /)
> lon2d(nf,0:203,:) = (/ a->Longitude /)
> else
> data(nf,:,:) = (/ x /)
> lat2d(nf,:,:) = (/ a->Latitude /)
> lon2d(nf,:,:) = (/ a->Longitude /)
> end if
> end if
> end do
>
> The script runs all the way through but I'm still finding a couple strange
> things.
> #1: The NetCDF file generated only shows ranges for Longitude :
> 85.1406-172.263 when using ncview, though there are values in there for the
> entire range...same as latitude, but since I know there are values there
> I'm not too concerned, just find it odd.
> #2: When I do plot there still looks like there is missing data, I've
> attached both again. Maybe it's removing duplicate values?
>
>
>
> -----Original Message-----
> From: Dennis Shea [mailto:shea at ucar.edu]
> Sent: Wednesday, May 27, 2015 12:43 PM
> To: BLIUJUS, STEVEN D CTR USAF AFMC AFLCMC/HBAW-OL
> Cc: ncl-talk at ucar.edu
> Subject: Re: [ncl-talk] Problems converting hdf to netcdf using NCL
>
> There is no reason for stopping after 34 iterations.
>
> ===
> [1] When I wrote "Each iteration thru the loop .... print the current
> 'time' value."
>
>
> I meant ...
>
> do nf=0,nfili-1 ; loop over files
> a = addfile(diri+fili(nf),"r")
> time(nf) = (/ parse_MODIS( fili(nf) ) /)
> print("nf="+nf+" time(nf)="+time(nf)) ; print current value
> ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
>
> [2] Also ... after the loop is complete ...
> ^^^^^
>
>
> end do
>
> ; manually add meta data
>
> time!0 = "time"
> time at long_name = "current time as YYYYDDDHHMM"
> time&time = time
> printVarSummary(time)
>
>
> print(time) ; print all values
>
>
>
> On Wed, May 27, 2015 at 10:08 AM, BLIUJUS, STEVEN D CTR USAF AFMC
> AFLCMC/HBAW-OL <steven.bliujus.3.ctr at us.af.mil> wrote:
>
>
> Dennis,
>
> Is there a default setting that loops kill after? The loop stops
> after 34 iterations, thus leaving out the other 226 files, which explains
> why there is no data anywhere else on the globe. There are so many files
> because we are doing a global analysis, so these are swaths through the
> entire day.
>
> As far as the plots.ncl. These values do match the values given by
> modis.ncl as far as lat/lon/aod. So the real problem is the loop quitting
> early. I attached the log file I created.
>
> Steve
>
> -----Original Message-----
> From: Dennis Shea [mailto:shea at ucar.edu]
>
> Sent: Wednesday, May 27, 2015 10:07 AM
> To: BLIUJUS, STEVEN D CTR USAF AFMC AFLCMC/HBAW-OL
> Cc: ncl-talk at ucar.edu
> Subject: Re: [ncl-talk] Problems converting hdf to netcdf using NCL
>
> I think you have to do some debugging.
>
> ===
>
>
> Look at "time": time: [20143050010..2.426116152246728e-315]
>
> The last time is obviously wrong! How did this happen?
>
> It did not happen with the 5 sample files ... why with all the
> files?
>
>
> Each iteration thru the loop .... print the current 'time' value.
>
>
> Also ... after the loop is complete ...
>
>
> print(time) ; print all values
>
>
> This issue must be cleared up!
>
>
> If 'time' is messed up ... what about other variables?
>
>
> -------
>
>
> In the code that reads/plots the netCDF ... did you do any
> debugging??
>
> diri = "./"
> fili = "MODIS_L2.Swath.nc"
> f = addfile(diri+fili, "r")
>
>
> print(f)
>
>
> AOT = f->AOT
>
> printVarSummary(AOT)
>
> printMinMax(AOT,0)
>
>
> lat2d = f->LAT2D
>
> printVarSummary(LAT2D)
>
> printMinMax(lat2d,0)
>
>
> lon2d = f->LON2D
>
> printVarSummary(lon2d)
> printMinMax(lon2d,0)
>
>
> Do these not match the original created values
>
>
> --------
>
> Minor: the dimension name "Cell_ACross_Swath" should be
> "Cell_Across_Swath"
>
>
> On Wed, May 27, 2015 at 8:34 AM, BLIUJUS, STEVEN D CTR USAF AFMC
> AFLCMC/HBAW-OL <steven.bliujus.3.ctr at us.af.mil> wrote:
>
>
> Dennis the script works however for some reason the data
> is not getting put into the NetCDF file correctly. I have ~120 files I'm
> going through, and when the script prints the dimensions, the following are
> stated:
>
> Variable: time
> Type: double
> Total Size: 2080 bytes
> 260 values
> Number of Dimensions: 1
> Dimensions and sizes: [time | 260]
> Coordinates:
> time: [20143050010..2.426116152246728e-315]
> Number Of Attributes: 1
> long_name : current time as YYYYDDDHHMM
>
> Variable: data
> Type: float
> Total Size: 28641600 bytes
> 7160400 values
> Number of Dimensions: 3
> Dimensions and sizes: [time | 260] x [Cell_Along_Swath |
> 204] x [Cell_ACross_Swath | 135]
> Coordinates:
> time: [20143050010..2.426116152246728e-315]
> Number Of Attributes: 2
> long_name : AOT at 0.55 micron for both ocean and land
> _FillValue : -9999
> (0) data: min=-0.05 max=4.277
>
> Variable: lat2d
> Type: float
> Total Size: 28641600 bytes
> 7160400 values
> Number of Dimensions: 3
> Dimensions and sizes: [time | 260] x [Cell_Along_Swath |
> 204] x [Cell_ACross_Swath | 135]
> Coordinates:
> time: [20143050010..2.426116152246728e-315]
> Number Of Attributes: 3
> units : degrees_north
> long_name : Geodetic Latitude
> _FillValue : -999
> (0) lat2d: min=-89.894 max=77.2616
>
> Variable: lon2d
> Type: float
> Total Size: 28641600 bytes
> 7160400 values
> Number of Dimensions: 3
> Dimensions and sizes: [time | 260] x [Cell_Along_Swath |
> 204] x [Cell_ACross_Swath | 135]
> Coordinates:
> time: [20143050010..2.426116152246728e-315]
> Number Of Attributes: 3
> units : degrees_east
> long_name : Geodetic Longitude
> _FillValue : -999
> (0) lon2d: min=-179.999 max=179.977
> (0) -----
>
>
> The thing that stands out is when the NetCDF file is
> created, the Longitude in the file only goes from 85.1406 - 172.263 and not
> -179.999 - 179.977 like the varSummary states. There is an error that gets
> thrown during the run:
>
> fatal:Dimension size mismatch on subscript #1, left-hand
> and right-hand side dimensions do not match
> fatal:["Execute.c":8128]:Execute: Error occurred at or
> near line 71 in file modis.ncl
>
> I have attached the image created from your script,
> compared to the image I get from the other script. I also reattached your
> script as I modified some directory things and added the libraries which
> shifted some of the lines. It should be noted that the information present
> from your script does match the data location on the image I created from a
> different script, it's just missing everything else.
>
>
> Steve
>
> -----Original Message-----
> From: Dennis Shea [mailto:shea at ucar.edu]
>
> Sent: Tuesday, May 26, 2015 5:43 PM
> To: BLIUJUS, STEVEN D CTR USAF AFMC AFLCMC/HBAW-OL
> Cc: ncl-talk at ucar.edu
> Subject: Re: [ncl-talk] Problems converting hdf to netcdf
> using NCL
>
> I just looked at the 5 sample files you sent ....
>
> The dimension sizes are *not* the same size across all
> files
>
>
> The 1st 4files are:
>
>
>
> MOD04_L2.A2014305.0000.006.2015077213911.pscs_000500937810.hdf
> dimensions:
> Cell_Along_Swath = 203 <=====
> Cell_Across_Swath = 135
>
> The last file is:
>
> MOD04_L2.A2014305.0500.006.2015077215106.pscs_000500937810
> dimensions:
> Cell_Along_Swath = 204 <====
> Cell_Across_Swath= 135
>
>
> That can be 'programmed around'.
>
>
> See attached
>
>
>
>
> On Tue, May 26, 2015 at 6:34 AM, BLIUJUS, STEVEN D CTR
> USAF AFMC AFLCMC/HBAW-OL <steven.bliujus.3.ctr at us.af.mil> wrote:
>
>
> Dennis,
>
> #2 and #3: Cell_Along_Swath, and Cell_Across_Swath
> are exactly the same for every file as far as actual dimension size, the
> variables are not though.
> #4: I am trying to create a mosaic of all the
> files and that is why I am wanting to generate one NetCDF file. So when I
> run my script the image above is what comes out. This image is what I want
> to create in a NetCDF file (Globe.PNG).
> #5: I did this because I was following a script
> online that I found where I needed to give the main variable(Optical Depth)
> the latitude/longitude dimensions. It used (::5) which what I just
> followed. When I try plotting using just (:,:), the plot take an extremely
> long time to plot(though there is much more data prevalent)...If it is one
> file it doesn’t take long, but when you add them all in the mosaic form, it
> takes a much more substantial amount of time to generate the image.
> #6: I'm not sure exactly what you mean by this as
> there is not time field in the file.
> #7: See 5 attached MODIS files
> #8: I am trying to create a NetCDF file for a
> coworker to use for his verification software. He only reads in NetCDF
> files or text files (though he would prefer NetCDF files)
>
> I also attached the script (with all previous
> commented lines removed) that plots the files and also creates a text table
> of the data.
>
> Any help would be appreciated.
>
> Steve
>
>
>
> -----Original Message-----
> From: Dennis Shea [mailto:shea at ucar.edu]
>
> Sent: Saturday, May 23, 2015 2:30 PM
> To: BLIUJUS, STEVEN D CTR USAF AFMC AFLCMC/HBAW-OL
> Cc: ncl-talk at ucar.edu
> Subject: Re: [ncl-talk] Problems converting hdf to
> netcdf using NCL
>
> I have looked at this and there are *lots* of
> issues.
>
> "... taking several hdf files (200-400) and trying
> to generate one NetCDF file."
>
>
> So ... **Conceptually**, you want
>
>
> AOT(time, Cell_Along_Swath, Cell_Across_Swath)
>
>
> Thing is your geographic coordinates are
> 2-dimensional. see below
>
>
>
>
> [1]
>
> The files contain swath data. From your dump:
>
> path:
> MOD04_L2.A2014306.2330.006.2015077230656.pscs_000500937810.hdf
> file global attributes:
> HDFEOSVersion : HDFEOS_V2.17
> StructMetadata_0 : GROUP=SwathStructure
> GROUP=SWATH_1
> SwathName="mod04"
>
>
> [2] Are the dimension sizes
>
> (Cell_Along_Swath,Cell_Across_Swath) =>
> (203,135)
>
>
> exactly the same for *all* files?
>
>
> [3] If [2] is satisfied .... Are the *values*
> contained within the following variables
>
> Latitude ( Cell_Along_Swath,
> Cell_Across_Swath)
>
> Longitude ( Cell_Along_Swath,
> Cell_Across_Swath )
>
>
> *exactly* the same for *all* files? ie:
> homogeneous space!
>
>
> My (limited) experience with satellite swath
> data is that the answer is "no"
>
>
> [4] If [2] *and* [3] are not satisfied, I see no
> advantage of creating a single netCDF file.
>
> You would end up with 200-400 sets of
> dimension sizes *and* 200-400 variables
>
> associated 200-400two-dimensional spatial
> variables.
>
>
> [5] I have no idea why you are sub-sampling
> [decimating] the array by 5:
> (::5,::5) the resolution from (203,135) to
> ??(41,27)??. This seems drastic!
>
> Sub-sampling (array decimation) is commonly
> used used when the resolution is (say) (2500,4000) and the data are
> reasonably smooth.
>
>
>
> [6] In any case you should be using
>
> ListSetType(a,"join")
>
>
> to add a record dimension ... which you
> should fill in with 'time'
>
> (parse file name) or some sequential number.
>
>
> Please read
>
> https://www.ncl.ucar.edu/Document/Functions/Built-in/addfiles.shtml
>
> Example 3 illustrates what "join" will do.
>
>
>
> [7] I do not have your ?MODIS? data file. Based
> upon your dump, I created a script that creates a netCDF file for each HDF
> file and (I hope) a sample graphic.
>
>
> [8] I am not sure of ultimate use of the data but
> from an NCL perspective, there is no reason to convert to netCDF.
>
>
> ===========================
>
>
> Please upgrade to 6.3.0
>
> http://www.ncl.ucar.edu/Download/
>
> ============================
>
>
>
>
>
> On Fri, May 22, 2015 at 9:59 AM, BLIUJUS, STEVEN D
> CTR USAF AFMC AFLCMC/HBAW-OL <steven.bliujus.3.ctr at us.af.mil> wrote:
>
>
> I wanted to provide an example of what is
> in the hdf file. See attached.
>
> Steve
>
>
> -----Original Message-----
> From: ncl-talk-bounces at ucar.edu [mailto:
> ncl-talk-bounces at ucar.edu] On Behalf Of BLIUJUS, STEVEN D CTR USAF AFMC
> AFLCMC/HBAW-OL
> Sent: Friday, May 22, 2015 10:36 AM
> To: ncl-talk at ucar.edu
> Subject: [ncl-talk] Problems converting
> hdf to netcdf using NCL
>
> I am taking several hdf files (200-400)
> and trying to generate one NetCDF file. I am able to plot it on a grid, but
> when I try to write to a NetCDF file, I get the following error:
>
> fatal:Number of dimensions in parameter
> (2) of (filedimdef) is (2), (1) dimensions were expected
> fatal:["Execute.c":8128]:Execute: Error
> occurred at or near line 82 in file File.ncl
>
>
> The data I am creating is AOD with two
> dimensions (lat/lon). I am uncertain as to why it is expecting one
> dimension. Below is my script. FYI, the original files had to be
> manipulated in order to make AOD have the lat/lon dimensions. Any help
> would be greatly appreciated.
>
>
> load
> "/home/bliujuss/ncl/lib/ncarg/nclscripts/csm/gsn_code.ncl"
> load
> "/home/bliujuss/ncl/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
> load
> "/home/bliujuss/ncl/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
> load
> "/home/bliujuss/ncl/lib/ncarg/nclscripts/csm/shea_util.ncl"
> begin
>
> diri =
> "/home/bliujuss/Verification_data/Test/"
> fili = systemfunc("ls "+diri+"*.hdf")
> a = addfiles(fili,"r")
>
> ;ListSetType(a,"cat")
>
>
> var1=a[:]->Optical_Depth_Land_And_Ocean(:,:)
> Optical_Depth_Land_And_Ocean =
> short2flt(var1)
>
> var2=a[:]->Longitude(:,:)
> Longitude=var2
>
>
> var3 = a[:]->Latitude(:,:)
> Latitude=var3
>
> dims =
> dimsizes(Optical_Depth_Land_And_Ocean)
> nrows = dims(0) ;203
> npixels = dims(1) ;135
>
> lat2d = new((/nrows,npixels/),float)
> lon2d = new((/nrows,npixels/),float)
>
> do i=0,npixels-1
> lat2d(:,i) = Latitude(:,i)
> end do
>
> do i=0,nrows-1
> lon2d(i,:) = Longitude(i,:)
> end do
>
> Optical_Depth_Land_And_Ocean at lat2d =
> lat2d(::5,::5) Optical_Depth_Land_And_Ocean at lon2d = lon2d(::5,::5)
>
>
> wks = gsn_open_wks("x11","plot_aod")
> gsn_define_colormap(wks,"NCV_banded")
> setvalues NhlGetWorkspaceObjectId()
> "wsMaximumSize": 300000000
> end setvalues
>
> res = True
> ; plot mods desired
> res at cnLinesOn = False
> ; turn off contour lines
> res at cnFillOn = True
> ; color plot desired
> res at cnLineLabelsOn = False
> ; turn off contour lines
> res at gsnAddCyclic = False
> ; non-global data
> res at gsnSpreadColors = True
> ; use full range of colormap
> ;res at cnFillMode = "RasterFill"
> ; turn on raster mode
> res at lbLabelAutoStride = True
> ; nice spacing of lables
> res at gsnMaximize = True
> ; blow up plot as much as poss
> ;print(min(lat2d))
> ;print(max(lat2d))
> ;print(min(lon2d))
> ;print(max(lon2d))
>
> plot =
> gsn_csm_contour_map_ce(wks,Optical_Depth_Land_And_Ocean(::5,::5),res)
>
> print(max(Optical_Depth_Land_And_Ocean(::5,::5)))
>
> nlat = dimsizes(lat2d(::5,::5))
> nlon = dimsizes(lon2d(::5,::5))
>
> diro = "/home/bliujuss/Plots/"
> filo = "example.nc"
> system("rm -f" + diro + filo)
> fout = addfile (diro + filo, "c")
>
> setfileoption(fout,"DefineMode",True)
>
> dim_names = (/"lat", "lon"/)
> dim_sizes = (/nlat, nlon/)
> dimUnlim = (/False, False/)
>
> filedimdef(fout,dim_names,dim_sizes,dimUnlim)
>
> filevardef(fout, "lat",
> typeof(lat2d(::5,::5)),getvardims(lat2d(::5,::5)))
> filevardef(fout, "lon",
> typeof(lon2d(::5,::5)),getvardims(lon2d(::5,::5)))
> filevardef(fout, "AOT",
> typeof(Optical_Depth_Land_And_Ocean(::5,::5)),getvardims(Optical_Depth_Land_And_Ocean(::5,::5)))
>
> filevarattdef(fout,"AOT",
> Optical_Depth_Land_And_Ocean(::5,::5))
> filevarattdef(fout,"lat", lat2d(::5,::5))
> filevarattdef(fout,"lon", lon2d(::5,::5))
>
> setfileoption(fout,"DefineMode",False)
>
>
> fout->Optical_Depth_Land_And_Ocean(::5,::5) =
>
>
> fout->(/Optical_Depth_Land_And_Ocean(::5,::5)/)
>
> fout->lat2d(::5,::5) = (/lat2d(::5,::5)/)
> fout->lon2d(::5,::5) = (/lon2d(::5,::5)/)
>
> end
>
>
> Steven Bliujus
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options,
> unsubscribe:
>
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options,
> unsubscribe:
>
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
>
>
>
>
>
>
>
>
>
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150527/5ea66657/attachment.html
More information about the ncl-talk
mailing list