function is_particle_field, name, pnames
@common_blocks.inc
pnames = *parti.names
ind = where(name[0] eq pnames)
if ind[0] ne -1 then return, 1 else return,0
end

function list_hdf5_field_names, file_name, COUNT=count, plusone=plusone
if (N_elements(plusone) eq 0) then plusone=0
Common options
; looks for all datasets in an enzo packed hdf5 file and returns the
; list of names
; Tom Abel July 2007
;        August 2007: now also works for the h5f movie format. 
    list = ""
    ;command = "h5ls -Svr " + file_name + " | grep 'Dataset'"
    command = "h5ls -r " + '"'+file_name+'"' + " | grep -c particle"
    if verbose then print, 'spawn command: ', command
    spawn, command, lines
    if verbose then print, 'returned ', lines, ' lines of output.'
    if (fix(lines[0]) gt 0) then begin
       command = "h5ls -r " + '"'+file_name+'"' + " | grep -m 1 particle | cut -f2 -d/"
       spawn, command, lines
       grid_name = lines[0]
       command = "h5ls -r " + '"'+file_name+'"' + " | grep " + grid_name + " | grep Dataset"
    endif else begin
       command = "h5ls -r "+ file_name+ " | head -n 100 | grep Dataset"
    endelse
    if verbose then print, 'spawn command: ', command
    spawn, command, lines
    if verbose then  print, 'returned ', N_elements(lines), ' lines of output.'
    if N_elements(lines) eq 0 then print, command, 'Warning: returned zero output!'

    secondslash = strpos(lines[0], '/', 1)
    firstgrid = strmid(lines[0], 0, secondslash+1)
;    lastgrid = strmid(lines[n_elements(lines)-1], 0, secondslash+1)
    lines = lines[where(strmatch(lines, firstgrid+'*'))]
    lines = strmid(lines, secondslash+1)
    ec = strpos(lines, 'Dataset')-1
    for i=0,N_elements(lines)-1 do $
       list = [list, strtrim(strmid(lines[i], plusone, ec[i]),2)]
    list = list[1:*]
    if list[0] eq 'Dark_Matter_Density' then begin
       he = list[0]
       list[0] = list[1]
       list[1] = he
    endif
    if ARG_PRESENT(COUNT) then count=N_elements(list)
    print, 'read ', N_elements(list), ' data labels from ', file_name
;breaki
    return, list
end

function return_grids_with_particles_hack, file_name
;return grid groups that really have particles
; not only where the index file says it does ... 
 command = getenv('FIRST') + '/TOOLS/strippedh5ls.sh '+ file_name + ' | grep _position_x'
 spawn, command, lines
 nl = N_elements(lines)
 for i=0,nl-1 do begin
    dl = strpos(lines[i], '/p', /reverse_search)
    en = strmid(lines[i], 0, dl+1)
    if i eq 0 then ls = en else ls = [ls, en]
 endfor
return, ls
end


function list_movie_field_names, file_name, COUNT=count
  file_name = strmid(file_name, 0, strpos(file_name, ':::'))
  ind = where(file_test(file_name, /READ) ne 0)
  if ind[0] eq -1 or file_name eq '' then begin
     print, 'list_movie_field_names: ', file_name, ' not readable. Guessing Density ..'
     return, 'Density'
  endif
;  command = getenv('FIRST') + '/TOOLS/strippedh5ls.sh '+ file_name +
;  ' | head -n 100'
  if file_test(file_name+'.idx') then begin 
     header = { pi: 0.0, root_delta: 0.0, nFields: 0B, staggering: 0B}
     header_length = 4+4+1+1
     openr, 1, file_name+'.idx'
     ;; Read header
     readu, 1, header
     if (abs(header.pi - 3.14159) gt 1e-5) then begin
        print, header.pi, format='("PI = ", F0)'
        print, "Endian wrong in the movie index? Reopening but byte swaping. Wish me luck."
        close, 1
        openr, 1, file_name+'.idx', /swap_endian
        filesize = (FSTAT(1)).size
        readu, 1, header
     endif
     list = strarr(header.nFields)
     if (header.nFields gt 0) then begin
        for i=0,header.nfields-1 do begin
           trash = '1234567890123456789012345678901234567890123456789012345678901234'
           readu, 1, trash
           list[i] = trash
;           print, trash
        endfor
     endif
     close, 1
  end else begin
     command = 'h5ls -Sr '+ file_name + ' | head -n 100'
     spawn, command, lines
     Datal = where(strmatch(lines, '*Dataset*'))
     nl = N_elements(datal)
     list = ['']
     c = 0
                                ; read the first and consecutive data line entries 
     while ((c eq 0) or ((Datal[min([c,nl-1])]-datal[max([c-1,0])]) eq 1))  do begin 
        cpos = Datal[c]
        dl = strpos(lines[datal[c]], 'Dataset')
        bl = strmid(lines[datal[c]], 0, dl)
        dash = strpos(bl, '/', /reverse_search)
        fieldname = strtrim(strmid(bl, dash+1,dl-dash),2)
        list= [list, fieldname]
        c += 1
     endwhile
     list = list[1:*]
  endelse
  if ARG_PRESENT(COUNT) then count=N_elements(list)
  print, 'read ', N_elements(list), ' data labels from ', file_name
  return, list
end

function list_field_names, fname
@common_blocks.inc
  data_format = determine_data_format_from_file_name(fname)
if data_format eq 1 then list = *parti.names
if data_format eq 2 then list = *parti.names
if data_format eq 3 then list = list_hdf5_field_names(fname, plusone=1)
if data_format eq 4 then list = list_hdf_field_names(fname)
if data_format eq 5 then list = list_hdf5_field_names(fname)
if data_format eq 6 then list = list_movie_field_names(fname)
if data_format eq 7 then list = list_flash_field_names(fname)
return, list
end

PRO  read_all_grid_hdf5, gin, fields, all,  VERBOSE = verbose, $
                         STRIDE=stride, START=startin, COUNT=count
; grids is an array of the grid_structure specifying the grids one wants
; files is an array of strings with the file_names in which each grid
; is to be found
; fields is a String array with the Attribute strings of the datasets
; to be read
; an array of the struct readall is returned

  readall        = {readall,                         $
                    file_name: "",                   $
                    Np_total: 0L,                    $
                    dims: [0,0,0],                   $
                    li: [0,0,0],                     $  ; left index
                    ri: [0,0,0],                     $  ; right index
                    N_dfields: 0B,                   $
                    DATA_FIELDS: STRARR(100),        $
                    Data:     PTRARR(100)   }
  
  stride_defined= n_elements(stride) gt 0
  start_defined = n_elements(startin) gt 0
  count_defined = n_elements(count) gt 0
  
  files = gin.baryon_file
  fnpos = strpos(files[0] , ':::')
  if fnpos[0] ge 0 then begin
     mgident = strmid(files, fnpos+3)
     files = strmid(files, 0, fnpos)
  endif else mgident = '/Grid'+string(gin.num , format='(I8.8)')+'/'

; sort them by file name so we only open & close every file just once
  oind    = sort(files)
  files   = files[oind] 
  grids   = gin[oind]
  mgident = mgident[oind]

  if start_defined then start = startin[oind, *]
  if stride_defined then stride = stride[oind, *]
  if count_defined then count = count[oind, *]

  revind = sort(oind)           ; index to get it back in incoming order
  ind = uniq(files)
  ufiles = files[ind]
  nf = N_elements(ufiles)  ; number of unique files
  ng = N_elements(grids)
  all = replicate(readall, ng)  ; array of readall to return data for every grid
  ind = [0, ind+1]

  if keyword_set(verbose) then print, 'read_all_grid_hdf5: reading ', FIELDS, $
                                      ' fields from ', ng, ' grids in', nf,' files.'
  if files[0] eq '' then BEGIN
      print, 'read_all_grid_hdf5: the file "', file_name, '" does not exist!'
      print, ' returning nothing useful'
      RETURN
  ENDIF
  
  N_df_requested = N_ELEMENTS(fields)
  special = ['x', 'y', 'z','dx','dy','dz','volume','mass']

  intd = where(file_test(ufiles, /READ) ne 0)
  if intd[0] eq -1 then begin
     print, 'read_all_grid_hdf5: ', ufiles, ' no readable file found. Returning.'
     return
  endif

  data_format = determine_data_format_from_file_name(ufiles[0])
  particlemassname = 'particle mass'
  particlemassname = 'particle_mass'
  
  cnt = 0L
  for fi=0L,nf-1 DO BEGIN
     file_name = ufiles[fi]
                                ; ok crazy hack: in case this filename
                                ; is a directory use the next file
                                ; ... 
     tind = where(strmatch(ufiles, '*/') eq 0)
     if strmatch(file_name, '*/') and tind[0] ne -1 then file_name =  ufiles[tind[0]] 


; next three lines are the hack
     if data_format eq 6 and total(is_particle_field(fields)) ne 0 then begin
         ls = return_grids_with_particles_hack(file_name)
         print, N_elements(ls), ' actual grids with particles'
     endif else ls = ['']
     
     ngtf = ind[fi+1]-ind[fi]   ; number of grids to read from this files
     if ngtf le 0 then break
     if keyword_set(verbose) then print, 'reading ', ngtf, ' grids from', file_name
     fid = H5F_OPEN(file_name)  
; Open the dataset  
     for j=0L,ngtf-1 do begin    ; read all the data for this grid
        if data_format eq 5 then gident = mgident[ind[fi]+j]
        if data_format eq 6 then gident = mgident[ind[fi]+j]
        if data_format eq 3 then gident = ''

        cg = grids[ind[fi]+j]
        ra = cg.end_index-cg.start_index + 1
        dd = double(cg.right_edge-cg.left_edge)/ra
        d = ra
        this_stride= (stride_defined) ? reform(stride[cnt,*]) : (intarr(n_elements(d)) + 1)
        this_start = (start_defined) ? reform(start[cnt,*]) : intarr(n_elements(d))
        this_count = (count_defined) ? reform(count[cnt,*]) : d[*]

        FOR i=0,N_elements(fields)-1 DO BEGIN
           field = fields[i]
           if field eq 'mass' then field = 'Density'
           if field eq 'particle massdensity' then field = particlemassname

           ; hack next line
           if data_format eq 6 and total(gident eq ls) eq 0 and  total(is_particle_field(field)) gt 0 then break


           if (where(field eq special))[0] ne -1 then begin ; this is a special field
               if count_defined then ra = ((count[cnt,*]) > 1 ) else ra=cg.dim > 1
               in = INDGEN(product(ra), /l64)
               data = dblarr(ra[0],ra[1],ra[2])
;print,ra, count[cnt,*]
;               print, 'ra',ra ,'dim',cg.dim
               indi = array_indices(ra,in, /dimensions)
               case field of
                   'x': data[*] = cg.Left_edge[0] + (0.5D + DOUBLE(this_start[0]+indi[0,*]))*dd[0]
                   'y': data[*] = cg.Left_edge[1] + (0.5D + DOUBLE(this_start[1]+indi[1,*]))*dd[1]
                   'z': data[*] = cg.Left_edge[2] + (0.5D + DOUBLE(this_start[2]+indi[2,*]))*dd[2]
                   'dx': data[*] = dd[0]
                   'dy': data[*] = dd[1]
                   'dz': data[*] = dd[2]
                   'volume': data[*] = product(dd[where(cg.dim gt 1)])
                   ELSE: data[*] = 1.
               endcase
;debs= size(data)
;if debs[0] ne 3 then print, debs, ra
           endif else begin     ; not special, read from file
              if data_format ne 3 then begin
                 found = 0
                 hdfgid = H5G_OPEN(fid, gident)
                 ngroup = H5G_GET_NUM_OBJS(hdfgid)
                 H5G_CLOSE, hdfgid
                 kl = ['']
                 for gm=0,ngroup-1 do begin
                    mn = H5G_GET_MEMBER_NAME(fid, gident, gm)
                    if mn eq field then found = 1
                    kl = [kl, mn]
                 endfor
                 if not found then begin
                    if field eq 'particle mass' then field = 'particle_mass'
                    if (where(field eq kl))[0] ge 0 then found = 1
                 endif
              endif else found = 1
               data = fltarr(1,1,1)
               if found then begin
                   dsetid = H5D_OPEN(fid, gident+strtrim(field,2))  
                   dataspace_id = H5D_GET_SPACE(dsetid)  
; Retrieve the dimensions 
                   d = H5S_GET_SIMPLE_EXTENT_DIMS(dataspace_id)  

                   ratio = this_count/this_stride > 1
;                   if total(this_start lt d) ne 3 then continue
                   dims = N_elements(d)
                   if (total(d lt this_start+ratio) gt 0 && $
                       is_particle_field(field) ne 1) then begin
                      print, 'hdf5read: trouble!', this_start, ratio
                      continue
                   endif
                  if n_elements(d) ne 1 then begin
                       H5S_SELECT_HYPERSLAB, dataspace_id, this_start[0:dims-1], ratio[0:dims-1], STRIDE=this_stride[0:dims-1], /RESET  
                       memory_space_id = H5S_CREATE_SIMPLE(ratio)  
                       data = H5D_READ(dsetid, FILE_SPACE=dataspace_id, MEMORY_SPACE=memory_space_id)  
                       H5S_CLOSE, memory_space_id  
                   endif else begin
                       this_stride=this_stride[0]
                       data = H5D_READ(dsetid, FILE_SPACE=dataspace_id)  
                   ENDELSE
                   H5S_CLOSE, dataspace_id  
                   H5D_CLOSE, dsetid
               endif else print, 'warning requested field not found in ', gident, field
            ENDELSE
           
           all[cnt].file_name   = gident
           all[cnt].np_total    = N_elements(data)
           s = size(data)
           if (s[0] le 2) then s[3] = 1
           if (s[0] eq 1) then s[2] = 1
           all[cnt].li[*]       = this_start[*]
           all[cnt].ri[*]       = (s[1:3]+this_start-1)[*]
           all[cnt].dims        = all[cnt].ri - all[cnt].li + 1
           all[cnt].N_Dfields   = N_df_requested
           all[cnt].data_fields[i] = fields[i]
           if fields[i] eq 'particle massdensity' then data *= product(dd)
           if fields[i] eq 'mass' then data *= product(dd[where(dd gt 0.)])
           all[cnt].Data[i]     = ptr_new(data, /no_copy)
           
        ENDFOR                  ; loop over fields
        cnt += 1
     ENDFOR                     ; loop over grids in this file
     if check_cancel() then break

     H5F_CLOSE, fid
  ENDFOR                        ; loop over unique files
; return information in the order grids and files were specified
  all = all[revind]

  RETURN
END