source: TOOLS/CMIP6_FORCING/SCENARIOS/AER_TROP_EMISSIONS/REGRID/netcdfwrite.pro @ 4007

; 1.11 01/10/07
;
; Takes an IDL structure and writes to NetCDF much as netcdf_read
; but from the opposite direction. 
;
; Dimension variables are identified as the first vector in the structure 
; of a certain length, and linked to all following vectors and multi-D 
; arrays with dimensions the same size. 
;
; The approach is simple and powerful but with one big flaw - arrays
; with non-unique dimension sizes (e.g. fltarr[36,36,5] rather than
; fltarr[36,37,5]) cannot be written unless:
;
;   a) A quick and dirty solution that breaks the NetCDF dimension
;      links: the DIM_ALL keyword is used.
;
;   b) A tedious but exact solution: the dimensions are explicilty
;      specified with the DIMENSIONS keyword.
;
; It will be possible to actually define which are dimensions etc... 
; ...in a later version... perhaps by defining "d_latitude" or "d_pressure"
; as indicating a dimension? Or if this is not flexible enough
; just write your own netcdf file following the Chapter 6 of the IDL 
; Scientific Data Formats manual
;
; Parameters
; ----------
;
; I:    file     string containing filename and path to be written
; I:    field    IDL structure to be written to file, arranged as 
;                in the example below
;
; Keywords
; --------
;
; CLOBBER - overwrite any existing file (otherwise, we fail with error
;       status if there's an exisitng file there)
; 
; DIM_ALL - treat all vectors as dimensions. This is needed in order to
;       write arrays non-unique dimension sizes, e.g. fltarr[36,36,5]
;
; DIMENSIONS - specify the exact dimensions and links. See below for
;      an example.
;
; ATTRIBUTES - pass in a structure containing the attributes for each 
;      dimension and variable, in the order they are found in the data 
;      structure we are writing to file. The structure should contain a
;      strarr for each attribute. The attribute name when written to
;      file will be the variable name in the struture, in lowercase.
;
; Example
; -------
;
; To add "units" and "long_name" attributes to each variable 
; for a partial compatability with the NetCDF standard, e.g. when 
; writing the 3D array/structure ozone:
;
;      ozone = {lon:fltarr(96),lat:fltarr(73),pressure:fltarr(50),$
;         ozone:fltarr(96,73,50)}
;
;      attributes = {units:strarr(4),long_name:strarr(4)}
;      attributes.units = ['degrees_east','degrees_north','Pa','kg kg-1']
;      attributes.long_name = ['Longitude','Latitude','Pressure',$
;         'Ozone']
;
;      netcdfwrite, 'filename', ozone, attributes=attributes
;
; On the unix command line "ncdump -h filename" produces:
;
;      netcdf filename {
;      dimensions:
;         LON = 96 ;
;         LAT = 73 ;
;         PRESSURE = 50 ;
;      variables:
;         float LON(LON) ;
;                 LON:units = "degrees_east" ;
;                 LON:long_name = "Longitude" ;
;         float LAT(LAT) ;
;                 LAT:units = "degrees_north" ;
;                 LAT:long_name = "Latitude" ;
;         float PRESSURE(PRESSURE) ;
;                 PRESSURE:units = "Pa" ;
;                 PRESSURE:long_name = "Pressure" ;
;         float OZONE(PRESSURE, LAT, LON) ;
;                 OZONE:units = "kg kg-1" ;
;                 OZONE:long_name = "Ozone" ;
;      }
;
; Note that NCDUMP reverses the order of the dimensions in multi-
; dimensional arrays compared to IDL.
;
; To explicitly specify the dimensions and linkages (not needed in most
; cases including this example: it's only needed if the rule of matching 
; dimensions by length breaks down because there are two dimensions of the
; same length):
;
;      dimensions = {isdim:intarr(4), links:intarr(5,4)}
;      dimensions.isdim =  [1,1,1,0]  ; (1=dimension, 0=variable)
;      dimensions.links = [[-1,-1,-1,-1,-1],[-1,-1,-1,-1,-1],$
;          [-1,-1,-1,-1,-1],[0,1,2,-1,-1]]
;
;      netcdfwrite, 'filename', ozone, attributes=attributes, $
;        dimensions=dimensions
;
;
;
; AJG 22/5/2003
;

; ------------------------------------------------------------------
;
; Add a non-string variable to the NetCDF file
;
; INPUTS
;
; var_type            - IDL type code - see IDL help for SIZE function
; ncid, name, dim_ids - parameters passed to ncdf_vardef (see IDL help here too)
;
; RETURNS
;
; NetCDF variable ID
;
function add_var, var_type, ncid, name, dim_ids

case var_type of
  ; byte
  1: var_id = ncdf_vardef(ncid, name, dim_ids, /BYTE)

  ; int
  2: var_id = ncdf_vardef(ncid, name, dim_ids, /SHORT)

  ; long
  3: var_id = ncdf_vardef(ncid, name, dim_ids, /LONG)

  ; float
  4: var_id = ncdf_vardef(ncid, name, dim_ids, /FLOAT)

  ; double
  5: var_id = ncdf_vardef(ncid, name, dim_ids, /DOUBLE)

  ; catch all (will probably catch us out some-day...)
  else: var_id = ncdf_vardef(ncid, name, dim_ids, /FLOAT) 
endcase

return, var_id

end

; ------------------------------------------------------------------
;
; Main function
;

pro  netcdfwrite, file, field, clobber=clobber, dim_all=dim_all, $
        attributes=attributes, dimensions=dimensions

   on_error, 2 ; on error, return to caller (could be improved
               ; at a later date to catch IO errors and close the
               ; file correctly).

   if ~keyword_set(clobber) then clobber = 0

   tag_names = tag_names(field)
   nTags = N_TAGS(field)

   if keyword_set(attributes) then begin
     ; Attributes are being supplied. Check they match the
     ; data structure also supplied.
     att_names = tag_names(attributes)
     nAtts = n_tags(attributes)
     nInAtt = lonarr(nAtts)
     for i=0,nAtts-1 do nInAtt[i] = n_elements(attributes.(i))
     if ~array_equal(nInAtt,replicate(nTags,nAtts)) then message, $
       'Each attribute must have a member for each data structure element'
   endif
  
   ncid = ncdf_create(file, clobber = clobber)   

   ; Run through all data "tags" to identify "dimensions" and
   ; "variables". For this to work, the data structure must always 
   ; have dimensions before normal variables  
   struct_dims, field, dimension_sizes, is_string_max, var_type, $
     nDimensions, dim_num, iDimLink, dim_all=dim_all, dimensions=dimensions

   dim_ids = intarr(nTags)
   var_ids = intarr(nTags)
   string_id = 0 ; a unique ID for any hidden string dimensions

   ; initialise file, defining the dimensions and variables
   ncdf_control, ncid, /fill
   for iTag=0, nTags-1 do begin
   
      if dim_num[iTag] eq 0 then begin
      
         ; Each dimension will have an associated variable
         ; containing the values at each point
         dim_ids[iTag] = ncdf_dimdef(ncid, tag_names[iTag], $
            dimension_sizes[iTag])

         var_ids[iTag] = add_var(var_type[iTag], ncid, tag_names[iTag], $
           dim_ids[iTag])
            
      endif else begin      
         if is_string_max[iTag] eq 0 then begin

           var_ids[iTag] = add_var(var_type[iTag], ncid, tag_names[iTag], $
             dim_ids[reform(iDimlink[iTag,0:dim_num[iTag]-1])])
            
         endif else begin
         
            ; Deal with the string array special case

            ; * Note that a string variable will gain an extra, hidden
            ; dimension when we write it to file because of NetCDF 
            ; limitations.
            
            ; Create a new dimension for the array of characters which forms
            ; a string, with its length being that of the longest string
            ; in the original array of strings. Don't bother defining
            ; a variable for the dimension.
            string_id = string_id+1
            string_dim = ncdf_dimdef(ncid, $
               'STRING'+string(format='(i0)', string_id), $
               is_string_max[iTag])

            ; Combine the "IDL" dimensions with the extra hidden
            ; string dimension (which must be the first, "fastest varying")
            string_dim_ids = $
               [string_dim, dim_ids[reform(iDimlink[iTag,0:dim_num[iTag]-1])]]            
            var_ids[iTag] = ncdf_vardef(ncid, tag_names[iTag], $
               string_dim_ids, /CHAR)    
                     
         endelse
      endelse
      
      ; Write any associated atributes
      if keyword_set(attributes) then begin
        for iAtt = 0, nAtts-1 do begin
          ncdf_attput, ncid, var_ids[iTag], strlowcase(att_names[iAtt]), $
            attributes.(iAtt)[iTag] 
        endfor
      endif

   endfor
   ncdf_control, ncid, /endef

   ; write the data 
   for iTag=0, nTags-1 do begin
   
      ncdf_varput, ncid, var_ids[iTag], field.(iTag)
      
   endfor

   ; close the file
   ncdf_close, ncid
   
end