pro zapbadval,silent,invar,outvar,table0,table1,table2,table3,table4,table5

;***ZAPBADVAL replaces badvals with values extrapolated from nearby good 
;***values. This is done using the average of usually 3 points from each end
;***of the gap.
;***
;***Input      silent = 1 for quiet operation, 0 for printing of run messages
;               invar = a vector or an array of data with badvals=-999.
;***                    Data must be in columns.
;***Output     
;***           outvar = an array exactly like invar, but it has been repaired
;***                    for holes by extrapolating between the endpoints.
;***                    The endpoints are generally based on the average of
;***                    3 points at each end of the gap. If only one or two
;***                    points exist, the first point next to the gap is used
;***                    for the extrapolation. If a gap is at the very
;***                    beginning or very end of invar, usually 3 points to
;***                    the right or left, respectively, are averaged and
;***                    the entire gap is given the avg. value. If only one or
;***                    two points are available, the first point to the
;***                    right and left, respectively, are used.
;***           tableX = a table of the badvalue gaps in invar. The X refers
;***                    to the column of invar that the table corresponds to.
;***                    format:[blocknumber,leftindex,rightindex]
;***                    (if no gaps, table=-1) 
;***                    Currently, invar arrays with any number of columns 
;***                    will get repaired, but only six tables can be
;***                    returned - this can be remedied by adding tables to
;***                    the first line of the code.
;***
;***Written by BWBerger 9/26/98 (Carved out of calrunavg.pro)

;***Initiallize an outvar the same size as invar
outvar=fltarr(col(invar),row(invar))-999.

;***Loop through the program once for each column of invar
for ii=0,col(invar)-1 do begin

if silent eq 0 then print,'Working on column ',strtrim(string(ii),2)

;***Define an r vector containing all the indicies of bad points in invar
r=where(invar(ii,*) eq -999.,c)

if (c ne 0) and (c ne row(invar)) then begin

;***Create a table of all the bad blocks of points in r
if len(r) eq 1 then begin
   if silent eq 0 then print,'only one bad pt in series'
   table=[0,r(0),r(0)]
endif else begin
   table=lonarr(3,50000)
   i=1L
   done=0
   blockcount=1L
   blocksize=1L
   while (done ne 1) do begin
      if r(i) eq r(i-1)+1 then begin
            blocksize=blocksize+1
         if (i eq len(r)-1) then begin
            done=1
            ;add last blk to table
            table(*,blockcount-1)=[blockcount,r(i-blocksize+1),r(i)]
            table=table(*,0:blockcount-1)
         endif else begin
            i=i+1
         endelse
      endif else begin

         table(*,blockcount-1)=[blockcount,r(i-blocksize),r(i-1)]
         if i eq len(r)-1 then begin
            table(*,blockcount)=[blockcount+1,r(i),r(i)]
            table=table(*,0:blockcount)
            done=1
         endif else begin
            blockcount=blockcount+1
            blocksize=1L
            i=i+1
         endelse
      endelse    
   endwhile
endelse

;***Check first block for low and high end points to zapextrapolate
;***across the block
repairtype=0
;check low end of first block
if table(1,0) ge 3 then begin       ;at least 3 low end pts to avg for extrap
   sl=2
   dl=3
   if silent eq 0 then print,'First block - at least 3 low end pts'
endif
if (table(1,0) eq 1) or (table(1,0) eq 2) then begin   ;only 1 or 2 low end pts
   sl=1
   dl=1
   if silent eq 0 then print,'First block - only 1 or 2 low end pts'
endif
if table(1,0) eq 0 then begin    ;no low end pts
   sl=0
   dl=0
   repairtype=1
   if silent eq 0 then print,'First block - no low end pts'
endif
;check high end of first block
if row(table) gt 1 then begin     ;there's more than 1 block
   if silent eq 0 then print,'There is more than one block.'
   if table(1,1)-table(2,0) ge 4 then begin  ;enough high end pts
      sh=2
      dh=3
      if silent eq 0 then print,'First block - at least 3 high end pts'
   endif
   if (table(1,1)-table(2,0) eq 3) or (table(1,1)-table(2,0) eq 2) then begin
      sh=1                                                 ; only 1 or 2 pts
      dh=1
      if silent eq 0 then print,'First block - only 1 or 2 high end pts'
   endif
endif else begin  ;First block is also the last block
   if table(2,0) le row(invar)-1-3 then begin    ;enough high end pts
      sh=2
      dh=3
      if silent eq 0 then print,'First & Last block - at least 3 high end pts'
   endif
   if (table(2,0) eq row(invar)-1-2) or $
                                  (table(2,0) eq row(invar)-1-1) then begin 
      sh=1                                  ;only 1 or 2 pts
      dh=1
      if silent eq 0 then print,'First & Last block - only 1 or 2 high end pts'
   endif
   if table(2,0) eq row(invar)-1 then begin      ;no high end pts
      sh=0
      dh=0
      repairtype=2
      if silent eq 0 then print,'First & Last block - no high end pts'
   endif
endelse
;repair invar
zapextrapo,invar(ii,*),table(1,0),table(2,0),sl,sh,dl,dh,repairtype,outvartmp
if row(table) eq 1 then begin
   goto, finish
endif

;***If there was more than one block, must now check last block.
repairtype=0
;check low end of last block
if table(1,row(table)-1)-table(2,row(table)-1-1) ge 4 then begin    
   sl=2                               ; enough low end pts
   dl=3
   if silent eq 0 then print,'Last block - at least 3 low end pts'
endif
if (table(1,row(table)-1)-table(2,row(table)-1-1) eq 2) $
          or (table(1,row(table)-1)-table(2,row(table)-1-1) eq 3) then begin
   sl=1                            ;only 1 or 2 low end pts
   dl=1
   if silent eq 0 then print,'Last block - only 1 or 2 low end pts'
endif
;check high end of last block
if table(2,row(table)-1) le row(invar)-1-3 then begin  ;enough high end pts
   sh=2
   dh=3
   if silent eq 0 then print,'Last block - at least 3 high end pts'
endif
if (table(2,row(table)-1) eq row(invar)-1-2)$
        or (table(2,row(table)-1) eq row(invar)-1-1) then begin
   sh=1                                      ;only 1 or 2 high end pts
   dh=1
   if silent eq 0 then print,'Last block - only 1 or 2 high end pts'
endif
if table(2,row(table)-1) eq row(invar)-1 then begin     ;no high end pts
   sh=0
   dh=0
   repairtype=2
   if silent eq 0 then print,'Last block - no high end pts'
endif
;repair last block
zapextrapo,outvartmp,table(1,len(table)-1),table(2,len(table)-1),sl,sh,dl,dh,$
        repairtype,outvartmp
if row(table) eq 2 then begin         ;only 2 blocks to fix
   if silent eq 0 then print,'Only two blocks'
   goto, finish
endif

;***If more than two blocks, check middle blocks
repairtype=0
for i=1,row(table)-2 do begin
   ;check low end of block
   if table(1,i)-table(2,i-1) ge 4 then begin   ; enough low end pts
      sl=2
      dl=3
      if silent eq 0 then print,'Middle block - at least 3 low end pts'
   endif
 if (table(1,i)-table(2,i-1) eq 3) or (table(1,i)-table(2,i-1) eq 2) then begin
      sl=1                       ;only 1 or 2 pts
      dl=1
      if silent eq 0 then print,'Middle block - only 1 or 2 low end pts'
   endif
   ;check high end of block
   if table(1,i+1)-table(2,i) ge 4  then begin ;enough high end pts
      sh=2
      dh=3
      if silent eq 0 then print,'Middle block - at least 3 high end pts'
   endif
 if (table(1,i+1)-table(2,i) eq 3) or (table(1,i+1)-table(2,i) eq 2) then begin
      sh=1                        ; only 1 or 2 pts
      dh=1
      if silent eq 0 then print,'Middle block - only 1 or 2 high end pts'
   endif
   ;repair middle block
   zapextrapo,outvartmp,table(1,i),table(2,i),sl,sh,dl,dh,repairtype,outvartmp
endfor


endif else begin
   outvartmp=invar(ii,*)
   if c eq 0 then begin
      table=-1
      if silent eq 0 then print,'No gaps.'
   endif else begin
      table=-2
      if silent eq 0 then print,'No good points in series.'
   endelse
endelse

finish:
outvar(ii,*)=outvartmp
cmd="table"+strtrim(string(ii),2)+"=table"
rr=execute(cmd)

endfor


return
end