;;;;;;;;;;;;;;;;;;;;
;;;   rall.pro   ;;;
;;;;;;;;;;;;;;;;;;;;

;;;
;;;  Author: nils/ adapted from rall.pro 20-nov-2003    
;;;  $Id$
;;;
;;;  Description:
;;;   Read data from all processors in a xy-plane (select plane by 
;;;   setting 'target_ipz' (default=0)) and combine them into one array
;;;   for each variable.
;;;   Overwrites nx, ny, etc, from start.pro, thus you need to run
;;;   start.pro after this if you want to continue working on data
;;;   from individual processors (e.g. if you want to run r.pro).
;
function param2
; Dummy to keep IDL from complaining. The real param2() routine will be
; compiled below
end
;
;  need to run start first: check whether this has been done
;
if (n_elements(started) le 0) then begin
  message, "You need to run start.pro first: use `.rnew start'"
endif
;
;  Read startup parameters
;
@readstartpars

;
;  read global sizes
;
nprocx=0L & nprocy=0L & nprocz=0L
close,1
openr,1,datatopdir+'/'+dimfile
readf,1,mx,my,mz,nvar,naux
readf,1,prec
readf,1,nghostx,nghosty,nghostz
readf,1,nprocx,nprocy,nprocz
close,1
;
ncpus = nprocx*nprocy*nprocz
;
;  read local sizes
;
datadir=datatopdir+'/proc0'
mxloc=0L & myloc=0L & mzloc=0L
;
close,1
openr,1,datadir+'/'+dimfile
readf,1,mxloc,myloc,mzloc
close,1


;
;  must set mz=mzloc such that varcontent defines the correct arrays
;  allow this value to be read in interactively if equal to default.
;
mz=mzloc
default,target_ipz,0
if target_ipz eq 0 then read,'give target_ipz> ',target_ipz


;
nxloc=mxloc-2*nghostx
nyloc=myloc-2*nghosty
nzloc=mzloc-2*nghostz
;
;  read data
;
default, dataopdir, 'data'
default, varfile, 'var.dat'
default, ldustvelocity, 0
default, ldustdensity, 0
;
x = fltarr(mx) & y = fltarr(my) & z = fltarr(mz)
xloc = fltarr(mxloc) & yloc = fltarr(myloc) & zloc = fltarr(mzloc)

;
;  Read data
;
varcontent=pc_varcontent()
totalvars=(size(varcontent))[1]

; Prepare for read
readstring=''
content=''
for i=0L,totalvars-1L do begin
  readstring = readstring + ',' + varcontent[i].idlvarloc
  content    = content + ', ' + varcontent[i].variable
  ; Initialise variable
  if (varcontent[i].variable eq 'UNKNOWN') then $
           message, 'Unknown variable at position ' + str(i)  $
                                    + ' needs declaring in varcontent.pro', /INFO   
  if (execute(varcontent[i].idlvar+'='+varcontent[i].idlinit,0) ne 1) then $
           message, 'Error initialising ' + varcontent[i].variable $
                                    +' - '+ varcontent[i].idlvar, /INFO
  if (execute(varcontent[i].idlvarloc+'='+varcontent[i].idlinitloc,0) ne 1) then $
           message, 'Error initialising ' + varcontent[i].variable $
                                    +' - '+ varcontent[i].idlvarloc, /INFO
;If it's a vector quantity skip the required number of elements
  i=i+varcontent[i].skip
end

content = strmid(content,2)
;
for i=0,ncpus-1 do begin        ; read data from individual files
  tag='proc'+str(i)
  datadir=datatopdir+'/proc'+strtrim(i,2)
  ; read processor position
  dummy=''
  ipx=0L &ipy=0L &ipz=0L
  close,1
  openr,1,datadir+'/'+dimfile
  readf,1, dummy
  readf,1, dummy
  readf,1, dummy
  readf,1, ipx,ipy,ipz
  ; read data
  if ((i eq 0) and (quiet le 2)) then begin
    print,'File contains: '+content
    print, FORMAT='(A,$)', "Reading: "
  endif

  if (ipz eq target_ipz) then begin
      close,1
      openr,1, datadir+'/'+varfile, /F77
      if (quiet le 2) then print, FORMAT='(A," ",$)', tag
      if (execute('readu,1'+readstring) ne 1) then $
        message, 'Error reading: ' + 'readu,1'+readstring
      ;
      ;  read deltay in case of shear
      ;
      if (lshear) then begin
          readu,1, t, xloc, yloc, zloc, dx, dy, dz, deltay
      end else begin
          readu,1, t, xloc, yloc, zloc
      end
      
      close,1
      ;
      ;  Don't overwrite ghost zones of processor to the left (and
      ;  accordingly in y and z direction makes a difference on the
      ;  diagonals)
      ;
      if (ipx eq 0) then begin
          i0x=ipx*nxloc & i1x=i0x+mxloc-1
          i0xloc=0 & i1xloc=mxloc-1
      endif else begin
          i0x=ipx*nxloc+nghostx & i1x=i0x+mxloc-1-nghostx
          i0xloc=nghostx & i1xloc=mxloc-1
      endelse
      ;
      if (ipy eq 0) then begin
          i0y=ipy*nyloc & i1y=i0y+myloc-1
          i0yloc=0 & i1yloc=myloc-1
      endif else begin
          i0y=ipy*nyloc+nghosty & i1y=i0y+myloc-1-nghosty
          i0yloc=nghosty & i1yloc=myloc-1
      endelse
      ;
      i0z=0 & i1z=i0z+mzloc-1
      i0zloc=0 & i1zloc=mzloc-1
      ;
      x[i0x:i1x] = xloc[i0xloc:i1xloc]
      y[i0y:i1y] = yloc[i0yloc:i1yloc]
      z[i0z:i1z] = zloc[i0zloc:i1zloc]

      for iv=0L,totalvars-1L do begin
          cmd =   varcontent[iv].idlvar $
            + "[i0x:i1x,i0y:i1y,i0z:i1z,*]=" $
            + varcontent[iv].idlvarloc $
            +"[i0xloc:i1xloc,i0yloc:i1yloc,i0zloc:i1zloc,*]"         
          if (execute(cmd) ne 1) then $
            message, 'Error combining data for ' + varcontent[iv].variable         
          ; For vector quantities skip the required number of elements
          iv=iv+varcontent[iv].skip
      endfor
  endif
endfor
if (quiet le 2) then print

;
xx = spread(x, [1,2], [my,mz])
yy = spread(y, [0,2], [mx,mz])
zz = spread(z, [0,1], [mx,my])
rr = sqrt(xx^2+yy^2+zz^2)

; reset datadir to more reasonable default
datadir=datatopdir+'/proc0'
;
;  reset boundary values and nx,ny,nz for (full) physical domain (not
;  sub-domain)
;
l1=3 & l2=mx-4
m1=3 & m2=my-4
n1=3 & n2=mz-4
;
nx=mx-2*nghostx
ny=my-2*nghosty
nz=mz-2*nghostz
;
;  fix z3=ztop which was local top in start.pro
;
if (lgravz) then begin
  ztop=z[n2] & z3=ztop
endif
;
;  Summarize data
;
@varcontent_stats
;
;  free memory
;
undefine, uu_loc
undefine, lnrho_loc
undefine, ss_loc
undefine, aa_loc
undefine, lncc_loc
;
read_all = 1                    ; marker for r.pro

end

; End of file