;
; Computes rbar for each region, for ring widths
;
; Get ring width timeseries and locations
;
ncid=ncdf_open('tree_rwid_nh.nc')
ncdf_diminq,ncid,'time',dummy,ntime
ncdf_diminq,ncid,'station',dummy,nstat
ncdf_varget,ncid,'year',x
ncdf_varget,ncid,'width',density
ncdf_attget,ncid,'width','missing_value',valmiss
ncdf_varget,ncid,'fraction',weight
ncdf_close,ncid
;
; Reduce length of it
;
help,density,x
keeplist=where(x ge 1601,nkeep)       ; use 1856 for recent period
if nkeep gt 0 then begin
  density=density(keeplist,*)
  x=x(keeplist)
endif
help,density,x
;
misslist=where(density eq valmiss,nmiss)
density(misslist)=!values.f_nan
;
; Get lists of trees in each region
;
restore,filename='reglists.idlsave'
;
; Repeat for each region
;
allrbar=fltarr(nreg)
for i = 0 , nreg-1 do begin
  ;
  ; Extract number of chronologies and the timeseries
  ;
  print,'Doing region ',regname(i),' with ntree=',ntree(i)
  nsize=ntree(i)
  inx=density(*,treelist(0:nsize-1,i))
  ;
  ; Define a correlation matrix
  ;
  grandr=fltarr(nsize,nsize)
  ;
  ; Do each tree against all others (including itself)
  ;
  for j = 0 , nsize-1 do begin
    print,j
    iny=reform(inx(*,j))
    rowr=manycorr(inx,iny,/nan,nval=nval)     ; compute r values
    grandr(*,j)=rowr(*)
  endfor
  ;
  ; Total it and average it, but only do i<>j.
  ;
  totr=0.
  totn=0.
  for j = 1 , nsize-1 do begin
    for k = 0 , j-1 do begin
      if finite(grandr(j,k)) then begin
        totr=totr+grandr(j,k)
        totn=totn+1.
        if grandr(j,k) eq 1. then message,'Ooops!'
      endif
    endfor
  endfor
  totr=totr/totn
  ;
  allrbar(i)=totr
  ;
endfor
;
;save,filename='rbar_rwid_recent.idlsave',nreg,regname,allrbar
save,filename='rbar_rwid_full.idlsave',nreg,regname,allrbar
;
end