[ncl-talk] Trouble converting dat to netcdf

Stanley Edwin sedwin at alaska.edu
Sat Mar 19 15:43:01 MDT 2016


​
 70200.dat.gz
<https://drive.google.com/a/alaska.edu/file/d/0B5e9G7If5LpKNnpnYXc0a3NfODg/view?usp=drive_web>
​
I could not get this dat file to convert properly to NetCDF, files from NWS

Thank You

Stanley G. Edwin
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160319/62990b2f/attachment.html 
-------------- next part --------------
;----------------------------------------------------------------------
; This script reads the 7 columns of data off a 70261.dat file and 
; writes 6 columns of the data and the time stamps to a NetCDF file.
;----------------------------------------------------------------------

;---Read DAT file as strings so we can parse them.

  file_prefix  = "70200"
  dat_filename = file_prefix + ".dat"
  lines        = asciiread(dat_filename,-1,"string")

;---Get index values of lines that contain a "#"
  sects = str_match_ind(lines,"#")
  nsect = dimsizes(sects)

;
; Loop through each subsection and store the # of rows for 
; each section in an array "nrows".
;
; Also parse the timestamp and save the date info to 
; individual arrays.
;
  maxLevel = 0
  nrows = new(nsect,integer)
  yyyy  = new(nsect,integer)
  mm    = new(nsect,integer)
  dd    = new(nsect,integer)
  hh    = new(nsect,integer)
  do i=0,nsect-1
    tmp     = str_split(lines(sects(i))," ")
   ;print("tmp(0) = " + tmp(0) + ", tmp(1) = " + tmp(1))
    yyyy(i) = toint(str_get_cols(tmp(0),6,9))
    mm(i)   = toint(str_get_cols(tmp(0),10,11))
    dd(i)   = toint(str_get_cols(tmp(0),12,13))
    hh(i)   = toint(str_get_cols(tmp(0),14,15))
    nrows(i) = toint(tmp(1))
  end do
  npts = sum(nrows)
  print("Total # npts = " + npts)
  maxLevel = max(nrows)
  print("MaxLevel = " + maxLevel)

;---Construct a time array
  mn = new(nsect,integer)
  sc = new(nsect,double)
  mn = 0
  sc = 0
  units  = "hours since " + min(yyyy) + "-01-01 00:00:00"
  time   = cd_inv_calendar(yyyy,mm,dd,hh,mn,sc,units, 0)

  printMinMax(time,0)

;
; Loop through each section again, and this time
; read the data into the data_array variable.
;
; Each row of data is formatted into lined-up columns.
; For example:
;
;21 95700B  135  -290B-9999    0    0
;21100100B  135  -289B-9999  120   10
;
; Note that it is possible for columns to run into each other.
;
; Column #1: 00-01
; Column #2: 02-07 (column 8 is a letter)
; Column #3: 09-13 (column 14 is a letter)
; Column #4: 15-19 (column 20 is a letter)
; Column #5: 21-25
; Column #6: 26-30
; Column #7: 31-35
;

;---Create empty arrays to hold all columns of data.
  p    = new((/nsect, maxLevel/),float)
  h    = new((/nsect, maxLevel/),float)
  t    = new((/nsect, maxLevel/),float)
  dp   = new((/nsect, maxLevel/),float)
  wdir = new((/nsect, maxLevel/),float)
  wspd = new((/nsect, maxLevel/),float)

;
; Loop through each section of rows and parse the 
; seven columns of data.
;
  n = 0                 ; counter
  do i=0,nsect-1
    ibeg = sects(i)+1            ; sect(i) is the "#xxxxx" line
    iend = sects(i)+nrows(i)
    p(i, 0:nrows(i)-1)    = tofloat(str_get_cols(lines(ibeg:iend),  2,  7))
    h(i, 0:nrows(i)-1)    = tofloat(str_get_cols(lines(ibeg:iend),  9, 13))
    t(i, 0:nrows(i)-1)    = tofloat(str_get_cols(lines(ibeg:iend), 15, 19))
    dp(i, 0:nrows(i)-1)   = tofloat(str_get_cols(lines(ibeg:iend), 21, 25))
    wdir(i, 0:nrows(i)-1) = tofloat(str_get_cols(lines(ibeg:iend), 26, 30))
    wspd(i, 0:nrows(i)-1) = tofloat(str_get_cols(lines(ibeg:iend), 31, 35))
  end do
  
;---Make sure -9999 and -8888 values are flagged as missing.
  p    = where(p.eq.-9999.or.p.eq.-8888,p at _FillValue,p)
  h    = where(h.eq.-9999.or.h.eq.-8888,h at _FillValue,h)
  t    = where(t.eq.-9999.or.t.eq.-8888,t at _FillValue,t)
  dp   = where(dp.eq.-9999.or.dp.eq.-8888,dp at _FillValue,dp)
  wdir = where(wdir.eq.-9999.or.wdir.eq.-8888,wdir at _FillValue,wdir)
  wspd = where(wspd.eq.-9999.or.wspd.eq.-8888,wspd at _FillValue,wspd)

 ;print(nrows(0))
 ;print(p(0, :))
 ;print(h(0, :))
 ;print(t(0, :))
 ;print(dp(0, :))
 ;print(wdir(0, :))
 ;print(wspd(0, :))

;---Attach some metadata
  time!0    = "time"
  time&time = time          ; This makes "time" a coordinate array
  printVarSummary(time)

  levs = ispan(0, maxLevel-1, 1)
  levs!0    = "levs"
  levs&levs = levs          ; This makes "levs" a coordinate array
  printVarSummary(levs)
 ;print(levs)

  nrows!0         = "time"
  nrows&time      = time
  nrows at long_name = "Available Levels"
  nrows at units     = "None"

  p!0         = "time"
  p!1         = "levs"
  p at long_name = "pressure"
  p at units     = "Pa"

  h!0         = "time"
  h!1         = "levs"
  h at long_name = "height"
  h at units     = "m"

  t!0         = "time"
  t!1         = "levs"
  t at long_name = "temperature"
  t at units     = "degC"

  dp!0         = "time"
  dp!1         = "levs"
  dp at long_name = "dew point"
  dp at units     = "degC"

  wdir!0         = "time"
  wdir!1         = "levs"
  wdir at long_name = "wind direction"
  wdir at units     = "deg"

  wspd!0         = "time"
  wspd!1         = "levs"
  wspd at long_name = "wind speed"
  wspd at units     = "m/s"

;---Print out min/max of each column
  printMinMax(p,0)
  printMinMax(h,0)
  printMinMax(t,0)
  printMinMax(dp,0)
  printMinMax(wdir,0)
  printMinMax(wspd,0)

;----------------------------------------------------------------------
; Section to write data to a NetCDF file
;----------------------------------------------------------------------

;---Make sure the file doesn't already exist.
  cdf_filename = "NWS"+file_prefix + ".nc"
  system("/bin/rm -f " + cdf_filename)
  fout = addfile(cdf_filename, "c")

;---Write global attributes to the file (optional)
  fout at title         = "NCL Simple Approach to netCDF Creation"
  fout at Conventions   = "None"
  fout at creation_date = systemfunc ("date")

  filedimdef(fout, "time", -1, True)

;---Write data values (and metadata) to NetCDF file.
  fout->p    = p
  fout->h    = h
  fout->t    = t
  fout->dp   = dp
  fout->wdir = wdir
  fout->wspd = wspd
  fout->time = time
  fout->avlevs = nrows



More information about the ncl-talk mailing list