pro drip
;**************************************************************************
;This program will read all ECOM data directly in and process it for
;analysis by other IDL programs...
;  need to .run drip.pro to compile subroutine
;  need ibmtoieee2.so (which includes ibmtoieee2.c and error.c)
;**************************************************************************
;file1 = 'tapename'
;file2 = antap + '.log'
;file3 = caldta.dta
;file4 = gratfile
;file6 = antap + '.asp'

maxlns = 10
maxbns = 340
maxizbns = 36
nsec = 60
numlns = 0
numbns = 0
maxrgns = 4

w = intarr(7)
wstuff = intarr(7, maxlns)
wline = fltarr(maxlns)
wl1 = intarr(maxlns)
wh1 = intarr(maxlns)
wl = intarr(maxlns)
wh = intarr(maxlns)
wl2 = intarr(maxlns)
wh2 = intarr(maxlns)

sb1 = fltarr(maxlns)
sb = fltarr(maxlns)
sb2 = fltarr(maxlns)
dn1 = fltarr(maxlns)
dn2 = fltarr(maxlns)

s = fltarr(128)
wav = fltarr(128)
pbk1 = fltarr(maxlns,maxbns,36,2)
plin = fltarr(maxlns,maxbns,36,2)
pbk2 = fltarr(maxlns,maxbns,36,2)
pcnt = fltarr(maxlns,maxbns,36,2)

pttime = fltarr(maxbns)
ptglon = fltarr(maxbns,2)
ptglat = fltarr(maxbns,2)
ptblat = fltarr(maxbns,2)
ptszen = fltarr(maxbns,2)
ptzen = fltarr(maxbns,36,2)
ptazi = fltarr(maxbns,36,2)
ptsaz = fltarr(maxbns,36,2)

sqglon = fltarr(maxbns,2)
sqglat = fltarr(maxbns,2)
sqblat = fltarr(maxbns,2)
sqszen = fltarr(maxbns,2)
sqzen = fltarr(maxbns,36,2)
sqazi = fltarr(maxbns,36,2)
sqsaz = fltarr(maxbns,36,2)

sdglon = fltarr(maxbns,2)
sdglat = fltarr(maxbns,2)
sdblat = fltarr(maxbns,2)
sdszen = fltarr(maxbns,2)
sdzen = fltarr(maxbns,36,2)
sdazi = fltarr(maxbns,36,2)
sdsaz = fltarr(maxbns,36,2)

nup = intarr(maxrgns)
nsd = intarr(maxrgns)
ndn = intarr(maxrgns)
sc =  fltarr(128)
up = fltarr(maxrgns,128)
sd = fltarr(maxrgns,128)
dn = fltarr(maxrgns,128)
rtime = fltarr(maxrgns,2)
rglon = fltarr(maxrgns,2)
rglat = fltarr(maxrgns,2)
rblat = fltarr(maxrgns,2)
rszen = fltarr(maxrgns,2)

integrating = intarr(maxrgns)

;**************************************************************************
; subroutine getcor and block data
;**************************************************************************
forward_function cnvtint
corfactor=fltarr(128)
corfactor(29:127)=1.0
corfactor(9:28)=[0.2532,0.3817,0.4444,0.5405,0.5882,0.6452,0.6897,$
	0.7143,0.7576,0.7692,0.8000,0.8130,0.8333,0.8696,0.8696,$
	0.8929,0.9091,0.9259,0.9259,0.9434]

;**************************************************************************
; subroutine prhead (10-20) - read in header information 
;**************************************************************************
tape='tapename'
nsec = 0

read, 'ENTER THE TIME BIN IN SECONDS (0 = 60 - DEFAULT): ', nsec
print, ' '
if ((n_elements(nsec) ne 1) and (nsec gt 0)) then nsec = 60

read, 'ENTER THE ZENITH ANGLE BIN WIDTH (MIN WIDTH = 5): ', izbnw
print, ' '

nzbns = fix(180/izbnw)
if (nzbns gt maxizbns) then begin
   print, '***** NOT ENOUGH ZENITH ANGLE BINS *****'
   print, 'EDIT PROGRAM TO SET MAXBNS TO AT LEAST ', nzbns, $
          ' AND RECOMPILE'
   print, 'ABORTED'
   stop
endif

read,'ENTER THE NAME OF THE TAPE FILE TO BE PROCESSED: ',tape
openr,1,tape
nstuff=lonarr(4)
junk=bytarr(4)
readu,1,junk
readu,1,nstuff       ; # headers, # scans, # events, # em?
nstuff=swap_endian(nstuff)
junk=bytarr(8)
readu,1,junk

temp='x This will read in all the initial info : do not shorten xxx' 
;61 characters
readu,1,temp

ecom=strmid(temp,0,12)
if (ecom ne '78-1ECOM-721') then begin
   print,' ** THIS IS NOT AN ECOM TAPE! **'
   return
endif

antap=strlowcase(tape)
rev=strmid(temp,23,4)
year=strmid(temp,27,4)
day=strmid(temp,31,4)
tstart=long(strmid(temp,35,6))
tend=long(strmid(temp,41,6))
rate=strmid(temp,47,8)
scan=strmid(temp,55,6)
junk=bytarr(33)
readu,1,junk

openw,2,antap+'.log'

printf,2, 'ANALOG TAPE NUMBER ',antap
printf,2, 'REV. # ', rev, '  YEAR ', year, '  DAY ', day
printf,2, 'START TIME(SEC)', tstart
printf,2, '  END TIME(SEC)', tend
printf,2, 'TOTAL SCAN EVENTS', scan
printf,2, ' '
printf,2,'OPTIONAL RECORDS:'

junk='xxxxxxxxxxxxxxxxxxxxxxxxxxxxxx'
opt=junk+junk+junk
junk=bytarr(4)

for i=0,nstuff(0)-2 do begin
   readu,1,junk
   readu,1,opt
   printf,2,opt
   readu,1,junk
endfor

printf,2,'END OF OPTIONAL RECORD.'
printf,2

print, ' '
print,'THIS IS TAPE ',antap
print,'REV #: ',rev,'  DAY: ',day,'  YEAR: ',year
print,'TAPE STARTS AT: ',tstart,' SECONDS'
print,'TAPE ENDS AT:   ',tend,' SECONDS'
print, ' '
read,'ENTER START TIME FOR ANALYSIS (-1 FOR BEGINNING): ', tml
read,'ENTER END TIME FOR ANALYSIS (-1 FOR END): ', tmh
if (tml eq -1.0) then tml=tstart
if (tmh eq -1.0) then tmh=tend

print, ' '
print,'*** ANALYSIS STARTS AT: ',tml,' SEC ***'
printf,2,'ANALYSIS STARTS AT: ',tml,' SEC'
print,'*** ANALYSIS ENDS AT:   ',tmh,' SEC ***'
printf,2,'ANALYSIS ENDS AT:   ',tmh,' SEC'

print, ' '
read,'ENTER GRATING STEP (-1 FOR TAPE VALUE): ',ig
print, ' '
read, 'ENTER THE NUMBER OF REGIONS FOR SPECTRUM PLOTS (4 MAXIMUM):  ',numrgns
print, ' '

numbns=((tmh-tml)/nsec) +1
if (numbns gt maxbns) then begin
   print, '***** NOT ENOUGH TIME BINS *****'
   print, 'EDIT PROGRAM TO SET MAXBNS TO AT LEAST ', numbns, $
          ' AND RECOMPILE'
   print, 'ABORTED'
   stop
endif

if (numrgns gt 1) then begin
for irgn = 0, numrgns-1 do begin
   print, 'REGION ', irgn+1
   print, 'ENTER START TIME, END TIME'  
   read, starttm, endtm

   rtime(irgn, 0) = starttm
   rtime(irgn, 1) = endtm
   print, ' '
endfor
endif

junk=lonarr(1)
readu,1,junk

;**************************************************************************
; subroutine prscan (30) - reads scan record of the tape
;**************************************************************************
scanb=dblarr(500)
scane=dblarr(500)

printf,2,' '
printf,2,' ' 
printf,2,'	SCAN RECORDS'
printf,2,'         Time    Sync.    VTCW(MSB) VTCW(LSB)  VTCW'
check1=0		; = j in old program...indexes scanb
check2=0		; = k in old program...indexes scane
i=lindgen(120)*6
isccnt=0
cnt=0

for j=0,nstuff(1)-1 do begin
   iscan=intarr(722)
   readu,1,iscan
   iscan=swap_endian(iscan)
   iscan=iscan(0:719)

   time=cnvtint(iscan(i))*65536+cnvtint(iscan(i+1))+double((cnvtint($
      iscan(i+2))*65536+cnvtint(iscan(i+3)))/1000.)
  
   if min(time) eq 0 then begin
      time=time(where(time ne 0,cnt))
   endif

   isccnt=isccnt+cnt
   insync=(iscan(i+4) lt 0)
   ivtcwm=insync*32768+iscan(i+4)
   vtcw=(ivtcwm)*65536+cnvtint(iscan(i+5))

   for k=0, n_elements(time)-1 do begin
      if (time(k) ne 0) then $
         printf,2,time(k),insync(k),ivtcwm(k),iscan(k*6+5),vtcw(k)
      if (k eq 0) or (time(k) gt 86400) then begin 
         goto,jump
      endif

      if (insync(k) eq 0) then begin
         if (check2 ge check1) then begin
	    scanb(check1)=time(k)-2
	    check1=check1+1
         endif
      endif
  
      if (insync(k) ne 1) or (check2 ge check1) then goto,jump

      scane(check2)=time(k)+2
      check2=check2+1

      jump: if (j gt 1000) and (k gt 1000) then begin
               print,'No. of scan points exceeds 1000'
            endif
   endfor
   printf,2,' '
   printf,2,'	TOTAL NUMBER OF SCAN EVENTS = ',isccnt
endfor
printf,2

junk=intarr(2)
readu,1,junk

;**************************************************************************
; subroutine prevnt (40) - calls event record of the tape
;**************************************************************************
icevnt=0 

printf,2,' '
printf,2,'	EVENT RECORDS'
printf,2,'          Time       Code  Status VTCW(MSB) VTCW(LSB)  VTCW'

for j=0,nstuff(2)-1 do begin
   ievent=intarr(840)
   readu,1,ievent
   ievent=swap_endian(ievent)
   i=indgen(120)*7
   time=cnvtint(ievent(i))*65536+cnvtint(ievent(i+1))+double((cnvtint($
      ievent(i+2))*65536+cnvtint(ievent(i+3)))/1000.)
   
   if min(time) eq 0 then time=time(where(time ne 0,cnt))

   icevnt=icevnt+cnt
   ievnt=ievent(i+6)/256
   istat=ievent(i+6) mod 256
   vtcw=cnvtint(ievent(i+4))*65536+cnvtint(ievent(i+5))
 
   for k=0,n_elements(time)-1 do begin
      printf,2,time(k),ievnt(k),istat(k),ievent(k*7+4),ievent(k*7+5),vtcw(k)  
   endfor
endfor

   itime = time(*)

printf,2,' ' 
printf,2,'	TOTAL NUMBER OF EVENTS = ',icevnt
printf,2

;**************************************************************************
; skip em data (50)
;**************************************************************************

junk='xxxxxxxxxxxx'
readu,1,junk
junk='xxxxxxxx'
for i=0,nstuff(3)-2 do begin
   readu,1,junk
   temp=opt  
   readu,1,temp
endfor

temp='xxxxxxxxxxxxxxxxxxxx'
temp=temp+temp+temp+temp+'xxxxxx'
readu,1,temp
junk=lonarr(2)
readu,1,junk

;**************************************************************************
; subroutine readdta (60) - reads the aspect part of the data
; subroutine readispec (80) - reads the spectra part of the data 
;**************************************************************************
itbn = 0
print, 'SEARCHING FOR TIME RANGE...'

data=fltarr(124)
ispec=intarr(1308)
i = 0
on_ioerror,nomore
eof1=0.0
good = 0
tgood = 0
bad = 0
tbad = 0 
index = 1
x = 0
zen = 0

scn = where(scanb gt 0, scancnt)
scanb = scanb(0:scancnt)
scanb(scancnt) = 0
scane = scane(0:scancnt)
scane(scancnt) = 0

while (not(eof(1))) and (not(eof1)) do begin
   number = 0
   break = 0
   ok = 0 
   readu,1,data

;***************************************************************************
; subroutine ibmtoieee2 - converts real number to ibm format to ieee format
;***************************************************************************
   for j=0,n_elements(data)-1 do begin
     data(j)=call_external('ibmtoieee2.so','fibmtoieee2',data(j),/d_value)
   end

;   data = shift(data,1) ;for boss0402
;   data(0) = data(24) ;for boss0402

   i=i+1
   junk=intarr(26) 
   readu,1,junk

   readu,1,ispec
   ispec=swap_endian(ispec)
   junk=fltarr(8)
   readu,1,junk
;   ispec = shift(ispec,2) ;for boss0402

;**************************************************************************
; skip data not in desired time range
;**************************************************************************
   if (data(0) lt tml) then begin
      zen=90.0-data(123)
      number = 4
   endif

   if (data(0) gt tmh) then begin
      eof1 = 1.0
      number = 3
   endif

   if (number ne 4) and (number ne 3) then break = 1

;**************************************************************************
; subroutine crspc (80)
;**************************************************************************
  if (break eq 1) then begin 
   for j = 4, 89 do ispec(j) = ((ispec(j) * 3) + 1) / 2
  
   for j = 90, 131 do ispec(j) = ((ispec(j) * 3) + 2) / 5

;**************************************************************************
; skip bad data (90)
;**************************************************************************
   for k = 0, scancnt do begin 
       if (data(0) gt scanb(k)) and (data(0) lt scane(k)) then begin
          number = 1
          k = scancnt+1
          goto, cont
      endif

      if (data(0) lt scane(k)) then begin
         goto, shift
      endif

      shift: dss1 = (fix(ispec(1304)/8)) mod 2
             dss2 = (fix(ispec(1305)/64)) mod 2
             if (dss1 eq dss2) then begin
                printf, 2,'>** ECOM DSS ERROR! **'
                printf, 2,'> ** ERROR IGNORED **'
                print,'** ECOM DSS ERROR! **'
                print,' ** ERROR IGNORED **'
                number = 1
                k = scancnt+1
                goto, cont
             endif
      cont:  
   endfor

   if (number ne 1) then break = 2   
  endif

;***************************************************************************
; grating step and resolution mode
;***************************************************************************
  if (break eq 2) then begin
   igo = -1
   iro = -1
   ibl1 = -1
   ibl2 = -1 
   yn = 'x'

   ig2 = fix(2*float(ispec(1300))/14.57143+1.07)
   ir = fix((fix(ispec(1305)/8)) mod 8)

   if (ig2 eq 32) then begin
       ig2 = 33
   end
   if (ig eq -1) then ig = ig2
   if (ir lt 2) then ir = 1
   
   l1 = (ig eq igo) and (ir eq iro)
   l2 = (ibl1 eq ispec(1304)) and (ibl2 eq ispec(1305))
   l3 = (igo eq -1) and (iro eq -1)

   if (l1 and l2) then number = 2
   
   if (x ne 1) then begin
   if (l3) then begin
      if (ig ne ig2) then begin
         print, ' '
         print, '** INPUT GRATING STEP AND INITIAL GRATING STEP **'
         print, '          ** FROM TAPE DO NOT MATCH! **'
         print, 'INITIAL STEP FROM TAPE = ', ig2
         read, 'USE TAPE VALUE? (Y/N) ', yn
         if (yn eq 'Y') or (yn eq 'y') then ig = ig2
      endif

      print, ' '
      printf, 2, '> THE GRATING-RESOLUTION SETTING IS:'
      printf, 2, format ='(6x, i2, 6x, i2)', ig, ir
      print, 'THE GRATING-RESOLUTION SETTING IS:'
      print, format ='(6x, i2, 6x, i2)', ig, ir
      print, ' '

   endif else begin $
      if (not l1) then begin
         printf, 2, '>**THE GRATING RESOLUTION CHANGED TO: ',ig2,' ',ir,' !?!**'
         print, '**THE GRATING RESOLUTION CHANGED TO: ',ig2,' ',ir,' !?!**'
                 
      endif
      if (not l2) then begin
         printf, 2, '>** ECOM DSS FLAGS CHANGED !?! **'
         print, '** ECOM DSS FLAGS CHANGED !?! **'
      endif
          
      printf, 2, '>** CHANGED IGNORED **'
      print, '** CHANGED IGNORED **'
   endelse
   endif
   
   igo = ig2
   iro = ir
   ibl1 = ispec(1304)
   ibl2 = ispec(1305)
   if (not l3) then number = 2
   if (x eq 1) then number = 2

   if (number ne 2) then break = 3
  endif
    
;***************************************************************************
; subroutine convrt - get wavelength vs. sensitivity table
;***************************************************************************
  if (break eq 3) then begin 

   x = 1
   openr, 3, 'caldta.dta' 

   if (ir eq 1) then begin
      if (ig ne 16) and (ig ne 19) and (ig ne 23) then begin
         print, 'caldta.dta does NOT have the G-R combination as found'
         print, 'on the tape: ', ig,',', ir, ' Enter the correct value'
         read, 'Enter 0 0 to abort', ig, ir
         if (ig eq 0) and (ir eq 0) then stop
      endif
      c = 1+(ig-15)/4
   endif

   if (ir eq 2) then begin
      if (ig ne 10) and (ig ne 13) and (ig ne 16) and (ig ne 19) then begin
         print, 'caldta.dta does NOT have the G-R combination as found'
         print, 'on the tape: ', ig,',', ir, ' Enter the correct value'
         read, 'Enter 0 0 to abort', ig, ir
         if (ig eq 0) and (ir eq 0) then stop
      endif
      c = 4+(ig-10)/3
   endif
         
   if (ir eq 3) then begin
      if (ig ne 1) and (ig ne 2) and (ig ne 16) and (ig ne 23) and $
         (ig ne 33) then begin
         print, 'caldta.dta does NOT have the G-R combination as found'
         print, 'on the tape: ', ig,',', ir, ' Enter the correct value'
         read, 'Enter 0 0 to abort', ig, ir
         if (ig eq 0) and (ir eq 0) then stop
      endif

      if (ig le 2) then c = 7+ig
      if (ig gt 2) then c = fix(10+((ig-16)/7))
   endif

   c = (c-1)*32
   close, 3
   free_lun, 3

;   on_ioerror, err
   openr, 3, 'caldta.dta'

   for j = 0, c-1 do begin
      readf, 3, junk 
   endfor

;   err: print, '***** ERROR READING CALDTA.DTA *****'
;   on_ioerror, nomore

   readf, 3, grat, res
   readf, 3, format = '(16(8e12.4/))', wav
   readf, 3, format = '(16(8e12.4/))', s

   if (ir eq 1) then begin
      for k = 0, 3 do begin
         wav(k) = wav(4)
         wav(123+k+1) = wav(123)
      endfor
      for k = 32, 39 do begin
         wav(k) = wav(31)
      endfor
   endif

   if (ir eq 2) then begin
      for k = 0, 9 do begin
         wav(k) = wav(10)
         wav(117+k+1) = wav(117)
      endfor
   endif

   if (ir eq 3) then begin
      for k = 0, 7 do begin
         wav(k) = wav(8)    
         wav(119+k+1) = wav(118)
      endfor
      wav(119) = wav(118)
   endif

   if (grat eq ig) and (res eq ir) then begin
      print, ' '
      print, 'G-R COMBINATION FOUND OK'
   endif else begin $
      print, 'Grating Step and Resolution Read From caldta.dta Do NOT'
      print, 'Match the Values on Tape - Position caldta.dta Properly'
      print, 'Values on Tape - ', ig, ' ', ir
      print, 'Those on caldta.dta - ', grat, ' ', res
   endelse

   close, 3
   free_lun, 3

;***************************************************************************
; subroutine fopen and function gratfile (105 and 107)
;***************************************************************************
   gratfile = strcompress('grat.'+string(fix(ig)), /remove_all)
   openr, 4, gratfile
   g = 0

   on_ioerror, next
   while (not(eof(4))) do begin
      readf, 4, w
      wstuff(*,g) = w
      g = g+1   
   endwhile

   next: on_ioerror, nomore
   close, 4
   free_lun, 4

   count = 0
   nl = where(wstuff ne 0, count)
   numlns = count/7

   printf, 2, format = '("> ", i2, " LINES AVAILABLE WITH THIS GRATING STEP")', $
              numlns

   for nw = 0, numlns-1 do begin
      wline(nw) = wstuff(0, nw)
      wl1(nw) = wstuff(1, nw)
      wh1(nw) = wstuff(2, nw)
      wl(nw) = wstuff(3, nw)
      wh(nw) = wstuff(4, nw)
      wl2(nw) = wstuff(5, nw)
      wh2(nw) = wstuff(6, nw)

      print, nw+1, wline(nw)
      printf, 2, wline(nw), wl1(nw), wh1(nw), wl(nw), wh(nw), wl2(nw), wh2(nw)
      printf, 2 , ' '
 
      for m = wl1(nw), wh1(nw) do sb1(nw) = sb1(nw) + s(m-1)
      for m = wl(nw), wh(nw) do sb(nw) = sb(nw) + s(m-1)
      for m = wl2(nw), wh2(nw) do sb2(nw) = sb2(nw) + s(m-1)

      sb1(nw) = sb1(nw) / float(wh1(nw) - wl1(nw) + 1)
      sb(nw) = sb(nw) / float(wh(nw) - wl(nw) + 1)
      sb2(nw) = sb2(nw) / float(wh2(nw) - wl2(nw) + 1)

      if (sb1(nw) ne 0) then dn1(nw) = float(wh(nw) - wl(nw) + 1) / $
                            (float(wh1(nw) - wl1(nw) + 1) * 2 * sb1(nw))
      if (sb2(nw) ne 0) then dn2(nw) = float(wh(nw) - wl(nw) + 1) / $
                            (float(wh2(nw) - wl2(nw) + 1) * 2 * sb2(nw))
   endfor

   print, ' '
   read, 'ENTER THE NUMBER OF THE WAVELENGTH TO BE PROCESSED:  ', nline
   print, ' '
   print, 'TIME RANGE FOUND. NOW READING DATA...'

   number = 2
  endif 

;***************************************************************************
; find satellite position, zenith angle, and time of observation for 
; each spectrum (110)
;***************************************************************************
  if (number eq 2) then begin 
   for irep = 0, 9 do begin
      ozen = zen
      zen = 90. - data((irep+1)*10+23)
      azi = data((irep+1)*10+22)
      time = data((irep+1)*10+14)
      glon = data(11)
      glat = data(10)
      blat = double(asin(data(9)/data(3))*(360/(2*3.14)))
      szen = 90. - data(17)
      sazi = data(16)
      phase = double(cos(zen)*cos(szen)+sin(zen)*sin(szen)*cos(azi-sazi))

      itbn = (fix(time - tml)/nsec) + 1

      if (itbn gt maxbns) then itbn = maxbns
      if (itbn lt 1) then itbn = 1
 
      if (zen ge ozen) then begin inc = 0
      endif else begin inc = 1
      endelse

;print, zen, inc
;stop

      izbn = (fix(zen)/izbnw) + 1

      if (izbn gt maxizbns) then izbn = maxizbns
      if (izbn lt 1) then izbn = 1
      if (izbn le 18) then f1200 = 0.083 else $
         f1200 = 0.083 + 0.002 * float(izbn - 18)

      for nw = 0, numlns-1 do begin
         for p = 0, 2 do begin
            if (p eq 0) then begin
               n1 = wl1(nw)
               n2 = wh1(nw)
            endif
            if (p eq 1) then begin
               n1 = wl(nw)
               n2 = wh(nw)
            endif
            if (p eq 2) then begin
               n1 = wl2(nw)
               n2 = wh2(nw)
            endif 
         
            sum = 0
            for n = n1, n2 do begin
               ind = irep * 128 + n + 4
               sum = sum + float(ispec(ind-1) * corfactor(n-1))
            endfor

            if (p eq 0) then fint1 = sum
            if (p eq 1) then fint = sum
            if (p eq 2) then fint2 = sum
         endfor
         
         if ((wline(nw)-1200.0) eq 0) then begin
            fint1 = fint1*f1200
            fint2 = fint2*f1200
         endif

         pcnt(nw, itbn-1, izbn-1,inc) = pcnt(nw, itbn-1, izbn-1,inc) + 1    
         pbk1(nw, itbn-1, izbn-1,inc) = pbk1(nw, itbn-1, izbn-1,inc) + fint1 
         plin(nw, itbn-1, izbn-1,inc) = plin(nw, itbn-1, izbn-1,inc) + fint 
         pbk2(nw, itbn-1, izbn-1,inc) = pbk2(nw, itbn-1, izbn-1,inc) + fint2 
      endfor      

      for irgn = 0, numrgns-1 do begin
         if ((time lt rtime(irgn,0)) or (time gt rtime(irgn,1))) then begin
            number = 5
         endif else begin        
            if (integrating(irgn) ne 0) then begin
               rglon(irgn, 0) = glon
               rglat(irgn, 0) = glat
               rblat(irgn, 0) = blat
               rszen(irgn, 0) = szen
               integrating(irgn) = 1
            endif
         endelse
     
         if (number ne 5) then begin
            j = 4 + 128 * irep
            for k = 0, 127 do begin
               sc(k) = ispec(k+j)
            endfor
         
            if ((zen ge 0) and (zen lt 50)) then begin
               nup(irgn) = nup(irgn) + 1
               for ichan = 0, 127 do begin
                  up(irgn,ichan) = up(irgn,ichan)+sc(ichan)*corfactor(ichan)
               endfor
            endif

            if ((zen ge 70) and (zen lt 110)) then begin
               nsd(irgn) = nsd(irgn) + 1
               for ichan = 0, 127 do begin
                  sd(irgn,ichan) = sd(irgn,ichan)+sc(ichan)*corfactor(ichan)
               endfor
            endif
            
            if (abs(phase) lt 0.2588) then begin
               nsd(irgn) = nsd(irgn) + 1
               for ichan = 0, 127 do begin
                  sd(irgn,ichan) = sd(irgn,ichan)+sc(ichan)*corfactor(ichan)
               endfor
            endif

            if ((zen ge 130) and (zen le 180)) then begin
               ndn(irgn) = ndn(irgn) + 1
               for ichan = 0, 127 do begin
                  dn(irgn,ichan) = dn(irgn,ichan)+sc(ichan)*corfactor(ichan)
               endfor
            endif
         endif
          
         if (number eq 5) then begin
	    if (integrating(irgn) eq 0) then begin
               rglon(irgn, 1) = glon
               rglat(irgn, 1) = glat
               rblat(irgn, 1) = blat
               rszen(irgn, 1) = szen
               integrating(irgn) = 0
            endif      
         number = 2    
         endif
      endfor

      ptzen(itbn-1, izbn-1, inc) = ptzen(itbn-1, izbn-1, inc) + zen 
      sqzen(itbn-1, izbn-1, inc) = sqzen(itbn-1, izbn-1, inc) + zen*zen 
      
      if (azi lt -90.) and (ptazi(itbn-1, izbn-1, inc) gt 0) then $
         azi = azi + 360.
      if (azi gt 90.) and (ptazi(itbn-1, izbn-1, inc) lt 0) then $
         azi = azi - 360.
      ptazi(itbn-1, izbn-1, inc) = ptazi(itbn-1, izbn-1, inc) + azi 
      sqazi(itbn-1, izbn-1, inc) = sqazi(itbn-1, izbn-1, inc) + azi*azi 

      if (sazi lt -90.) and (ptsaz(itbn-1, izbn-1, inc) gt 0) then $
         sazi = sazi + 360.
      if (sazi gt 90.) and (ptsaz(itbn-1, izbn-1, inc) lt 0) then $
         sazi = sazi - 360.
      ptsaz(itbn-1, izbn-1, inc) = ptsaz(itbn-1, izbn-1, inc) + sazi 
      sqsaz(itbn-1, izbn-1, inc) = sqsaz(itbn-1, izbn-1, inc) + sazi*sazi 

      if (glon lt -90.) and (ptglon(itbn-1, inc) gt 0) then $
         glon = glon + 360.
      if (glon gt 90.) and (ptglon(itbn-1, inc) lt 0) then $
         glon = glon - 360.
      ptglon(itbn-1, inc) = ptglon(itbn-1, inc) + glon 
      sqglon(itbn-1, inc) = sqglon(itbn-1, inc) + glon*glon 
      ptglat(itbn-1, inc) = ptglat(itbn-1, inc) + glat 
      sqglat(itbn-1, inc) = sqglat(itbn-1, inc) + glat*glat 
      ptblat(itbn-1, inc) = ptblat(itbn-1, inc) + blat
      sqblat(itbn-1, inc) = sqblat(itbn-1, inc) + blat*blat 
      ptszen(itbn-1, inc) = ptszen(itbn-1, inc) + szen 
      sqszen(itbn-1, inc) = sqszen(itbn-1, inc) + szen*szen 
   endfor

   ok = 1
   number = 1
endif

;***************************************************************************
; reports good set or bad set and approximately how much data has been
; read (130) 
;***************************************************************************
  if (number eq 1) then begin
   if (not ok) then begin
      if (good ne 0) then begin
         printf, 2, 'READ', good, ' GOOD DATA SETS'
         good = 0
      endif
      bad = bad + 1
      tbad = tbad + 1
   endif 

   if (ok) then begin
      if (bad ne 0) then begin
         printf, 2, 'SKIPPED', bad, ' BAD DATA SETS'
         bad = 0
      endif 
      good = good + 1 
      tgood = tgood + 1 
   endif 

   if (itbn ge (numbns*index)/20) then begin
      print, format = '(i3, "% OF DATA HAS BEEN PROCESSED")', index*5
      index = index + 1
   endif 

   number = 0
  endif 

  if (number eq 4) then number = 0

endwhile
nomore: on_ioerror,null

number = 3

;***************************************************************************
; finished reading data - finds averages over each time bin (140)
;***************************************************************************
if (number eq 3) then begin

if (good ne 0) then printf, 2, 'READ', good, ' GOOD DATA SETS'
if (bad ne 0) then printf, 2, 'SKIPPED', bad, ' BAD DATA SETS'           

for itbn = 0, numbns-1 do begin
   pttime(itbn) = tml + itbn * nsec + (nsec/2)
   for inc = 0, 1 do begin
      count = 0
      for izbn = 0, nzbns-1 do begin
         count = count + pcnt(nline-1,itbn, izbn, inc)
      endfor   
      if (count ne 0) then begin
         if (count ne 1) then begin
            fcnt = float(count)
            sdglon(itbn, inc) = sqrt(abs((fcnt*sqglon(itbn,inc)-$
	 	ptglon(itbn,inc)^2)/fcnt/(fcnt - 1.0)))
	    sdglat(itbn, inc) = sqrt(abs((fcnt*sqglat(itbn,inc)-$
		ptglat(itbn,inc)^2)/fcnt/(fcnt - 1.0)))
	    sdblat(itbn, inc) = sqrt(abs((fcnt*sqblat(itbn,inc)-$
		ptblat(itbn,inc)^2)/fcnt/(fcnt - 1.0)))
	    sdszen(itbn, inc) = sqrt(abs((fcnt*sqszen(itbn,inc)-$
		ptszen(itbn,inc)^2)/fcnt/(fcnt - 1.0)))
         endif else begin $
            sdglon(itbn, inc) = 0
            sdglat(itbn, inc) = 0
            sdblat(itbn, inc) = 0
            sdszen(itbn, inc) = 0
         endelse
         ptglon(itbn, inc) = ptglon(itbn, inc)/float(count)
         ptglat(itbn, inc) = ptglat(itbn, inc)/float(count)
         ptblat(itbn, inc) = ptblat(itbn, inc)/float(count)
         ptszen(itbn, inc) = ptszen(itbn, inc)/float(count)
      endif
   endfor
endfor  

for irgn = 0, numrgns-1 do begin
   if (integrating(irgn)) then begin
      rglon(irgn, 1) = glon
      rglat(irgn, 1) = glat
      rblat(irgn, 1) = blat
      rszen(irgn, 1) = szen
   endif   
endfor

;**************************************************************************
; print data to file5 (all) and file6 (time, glon, glat, blat, szen)
;**************************************************************************
openw, 6, antap + '.asp'
printf, 6, antap, rev, day, year

for inc = 0, 1 do begin
   for itbn = 0, numbns-1 do begin
      printf, 6, pttime(itbn), ptglon(itbn, inc), ptglat(itbn, inc), $
                 ptblat(itbn, inc), ptszen(itbn, inc)

      for izbn = 0, nzbns-1 do begin
         fcnt = pcnt(nline-1,itbn, izbn, inc)
         if (fcnt ge 2) then begin
            sdzen(itbn, izbn, inc) = sqrt(abs((fcnt*sqzen(itbn, izbn, inc) -$
                                     ptzen(itbn, izbn, inc)^2)/fcnt/(fcnt-1.)))
            sdazi(itbn, izbn, inc) = sqrt(abs((fcnt*sqazi(itbn, izbn, inc) -$
                                     ptazi(itbn, izbn, inc)^2)/fcnt/(fcnt-1.)))
            sdsaz(itbn, izbn, inc) = sqrt(abs((fcnt*sqsaz(itbn, izbn, inc) -$
                                     ptsaz(itbn, izbn, inc)^2)/fcnt/(fcnt-1.)))
         endif else begin $
            sdazi(itbn, izbn, inc) = 0
            sdsaz(itbn, izbn, inc) = 0
         endelse
    
         if (fcnt ne 0) then begin
            ptzen(itbn, izbn, inc) = ptzen(itbn, izbn, inc)/fcnt
            ptazi(itbn, izbn, inc) = ptazi(itbn, izbn, inc)/fcnt
            ptsaz(itbn, izbn, inc) = ptsaz(itbn, izbn, inc)/fcnt
         endif
      endfor
   endfor
endfor

for itbn = 0, numbns-1 do begin
   for izbn = 0, nzbns-1 do begin 
      printf, 6, pttime(itbn), ptazi(itbn, izbn, 0), ptazi(itbn, izbn, 1)
   endfor
endfor


tml = long(tml)
tmh = long(tmh)
numbns = long(numbns)
numlns = long(numlns)
numrgns = float(numrgns)
itime = float(itime)

ix = where((ievnt eq 1) and (itime ge tml) and (itime le tmh) , cn)
if (cn gt 0) then begin
   ievnt = ievnt(ix)
   istat = istat(ix)
   itime = itime(ix)
endif

ix = fix(ix)

save, filename = antap+'.one', $
      string(antap),string(rev), string(day), string(year), tml, tmh,$
      numbns, numlns, numrgns, ig, fix(ir), nline, $
      float(wline), fix(wl1), fix(wh1), fix(wl), fix(wh), $
      fix(wl2), fix(wh2), float(sb1), float(sb), float(sb2), $
      float(dn1), float(dn2), float(pttime), float(ptglon), float(ptglat), $
      float(ptblat), float(ptszen), float(sdglon), float(sdglat), float(sdblat),$
      float(sdszen), float(pcnt), float(pbk1), float(plin), float(pbk2), $
      float(ptzen), float(sdzen), float(ptazi), float(sdazi), float(ptsaz), $
      float(sdsaz), float(s), float(wav), fix(nup), fix(nsd), fix(ndn), $
      float(rtime), float(rglon), float(rglat), float(rblat), float(rszen), $
      float(up), float(sd), float(dn), float(izbnw), fix(nsec), fix(ievnt),$
      fix(istat), itime, ix

print, ' '
print, 'FINISHED WITH THIS TAPE'
print, format = '("DATA STORED IN FILE", a10, ".one")', antap

printf, 2, ' '
printf, 2, 'FINISHED WITH THIS TAPE'
printf, 2, format = '("TOTAL NUMBER OF BAD DATA SETS: ", i5)', tbad
printf, 2, format = '("TOTAL NUMBER OF GOOD DATA SETS: ", i5)', tgood
printf, 2, format = '("DATA STORED IN FILE", a10, ".one")', antap

close, 6
free_lun, 6
endif

close, 1
free_lun, 1
close, 2
free_lun, 2  

return
end

;**************************************************************************
function cnvtint,intin
return, temp=65536*(intin lt 0)+intin
end

;**************************************************************************
;		EXPLANATION OF ARRAYS AND VARIABLES
;		-----------------------------------
;
;azi - line of sight azimuth
;
;blat - magnetic latitude
;
;caldta.dta - file that contains wavelength and sensitivity table for 
;   various grating steps
;
;corfactor - correction factor table
;
;data - aspect information
;
;dn1 & dn2 - normalization factors which take into account that the line
;   and background region may differ in size
;
;glat - geographic latitude
;
;glon - geographic longitude
;
;gratfile - file that contains the wavelengths of interest that can be seen 
;   with the chosen grating step and the bin numbers for the line width and
;   background regions 
;
;ig - grating step 
;
;ir - resolution mode
;
;irep - current spectrum
;
;ispec - contains count values for each of the ten spectra
;
;itbn - time bin
;
;izbn - zenith angle bin
;
;nw - line number 
;
;pcnt - number of spectra in each bin
;
;pbk1 - contains counts in background range one
;plin - contains counts in a line
;pbk2 - contains counts in background range two
;
;pt prefix - sum of
;
;s - sensitivity
;
;saz - solar azimuth angle
;
;sb1 - average sensitivity for background range one
;sb - average sensitivity of the line
;sb2 - average sensitivity for background range two
;
;scanb & scane - time ranges during which data is known to be bad
;
;sd prefix - standard deviation of
;
;sq prefix - sum of squared 
;
;szen - solar zenith angle
;
;tapename - file that contains original data
;tapename.one - file that contains processed line and spectrum data
;tapename.log - file that contains header information, status, number of
;   good/bad data sets read, etc.
;
;time - time of observation
;
;tml - start time
;
;tmh - end time
;
;wav - wavelength
;
;wline - line wavelength
;wl1 - low channel limits for low wavelengths
;wh1 - high channel limits for low wavelengths
;wl - low channel limits for line
;wh - high channel limits for line
;wl2 - low channel limits for high wavelengths
;wh2 - high channel limits for high wavelengths
;
;zen - line of sight zenith angle (look direction)
;
;**************************************************************************