; Procedure to compare SLODAR and MASS-DIMM profiles
; Input: month and date of the evening
; Dec 2, 2004
;---------------------------------------------------------
pro mdseeing2, mm, dd

; slodar data selection
restore, 'slodar.idl'
; ---------slodar structure:
; prof = { UTyear: 0, UTmonth: 0, UTday: 0, UT: 0.0, $
;	  r0: 0.0, alt:FLTARR(9), Intcn2:FLTARR(9) }

; wrong (evening) dates in the CN2 files!
eve = where((slodat.UTmonth eq mm) and (slodat.UTday eq dd) )
neve = n_elements(eve) 
morn = -1
nmorn = 0

; Correct UT dates
;eve = where((slodat.UTmonth eq mm) and (slodat.UTday eq dd) and (slodat.UThour gt 12) )
;dd1 = dd+1   & mm1 = mm
;if (dd1 gt 30) then begin dd1=1 & mm1=mm+1 & endif
;morn = where((slodat.UTmonth eq mm1) and (slodat.UTday eq dd1) and (slodat.UThour lt 12) )

;nmorn =  n_elements(morn)
nslo = neve + nmorn

if (nslo eq 0) then return
prof = slodat[0]
out = replicate(prof,nslo)

out(0:neve-1) = slodat(eve)
;out(neve:nslo-1) = slodat(morn)

utslo = fltarr(nslo)
for i=0,nslo-1 do utslo[i] = out[i].UT

print, 'Total SLODAR profiles found: ', nslo
print, 'UT range: ',utslo[0], utslo[nslo-1]
nalt = 8
tmp = fltarr(nslo,nalt) 
for i=0,nslo-1 do tmp(i,*) = out[i].intcn2(0:nalt-1)

;--------MASS-DIMM data --------------------------------
 restore, 'massdimm.idl'
; mddat structure:
; prof = { UTyear: 0, UTmonth: 0, UTday: 0, UTmass:0.0, $
;	  fsee:0.0, Intcn2:FLTARR(6), UTdimm:0.0, see:0.0 }

;eve = where((mddat.UTmonth eq mm) and (mddat.UTday eq dd) and (mddat.UTmass gt 12) )
dd1 = dd+1   & mm1 = mm
if (dd1 gt 30) then begin dd1=1 & mm1=mm+1 & endif

morn = where((mddat.UTmonth eq mm1) and (mddat.UTday eq dd1) and (mddat.UTmass lt 12) )

;neve = n_elements(eve)
neve = 0
nmorn =  n_elements(morn)
nmd = neve + nmorn

prof = mddat[0]
md = replicate(prof,nmd)
; md(0:neve-1) = mddat(eve)
md(neve:nmd-1) = mddat(morn)

utmd = fltarr(nmd)
for i=0,nmd-1 do utmd[i] = md[i].UTmass


;----------------------------------------------------
; Data analysis
;-------------------------------------------------------

hmass = [0.01,0.5,1]*1000.   ; MASS layers: ground, 0.5, higher
nmass = 3
fact = alog10(2.)


; MASS triangular weights
;wt = fltarr(nmass,nalt)
;for i=1,nmass-1 do   wt(i,*) = (1. - abs(alog10(realt/hmass[i])/fact) ) > 0.
;for i=1,nalt-1 do if (realt[i] gt 16000) then wt(nMASS-1,i) = 1.
;wt(0,*) = 1. - total(wt(1:nmass-1,*),1)


utslo = out.ut

;------- spatial binning of SLODAR data
sloprof = fltarr(nslo,nmass)

  for i=0,nslo-1 do begin
   dh = out[i].alt[1] -  out[i].alt[0]
   if (dh lt 500.) then begin
   wt = fltarr(nmass,9)
   wt(1,*) = (1. - abs(alog10(out[i].alt/hmass[1])/fact) ) > 0.
  low = where(out[i].alt le hmass[1])
   wt(0,low )  = 1. - wt(1,low)
   wt(2,*)  = ( 1 - wt(1,*)-wt(0,*) ) > 0.
  for k=0,nmass-1 do  sloprof(i,k) = total(out[i].intcn2*wt(k,*))*0.01

   endif  
endfor

print, 'Spatial binning done!'
;stop

; temporal binning
utstart = 1.
utend = 9.
utstep = 0.2
nut = fix(utend - utstart)/utstep
rebinut = findgen(nut)*utstep + utstart

reslo = fltarr(nut,nmass)
remd = fltarr(nut,nmass)
numslo = intarr(nut) & nummd = intarr(nut)

for k=0,nut-1 do begin
  for i=0,nslo-1 do begin
  if ((abs(out[i].ut - rebinut[k]) le utstep*0.5) and (sloprof[i,0] gt 0.)) then begin
    numslo[k] += 1
    reslo(k,*) += sloprof(i,*)
  endif
  endfor

  for i=0,nmd-1 do begin
  if (abs(md[i].UTmass - rebinut[k]) le utstep*0.5)  then begin
    nummd[k] += 1
    remd(k,1) += md[k].intcn2[0]
    remd(k,2) += total(md[k].intcn2(1:5) )
    remd[k,0] += ( (md[i].see^(5./3.) - md[i].fsee^(5./3.) ) > 0) *6.8
  endif
  endfor

endfor

sel = where ((nummd gt 1) and (numslo gt 1))
if (sel[0] lt 0) then begin
  print, 'No useful data found!'
  return
endif

rebinut = rebinut(sel)
nummd = nummd(sel) & numslo = numslo(sel)
nut = n_elements(rebinut)
for k=0,nut-1 do begin
  remd(k,*) = remd(sel[k],*)/nummd[k]  
  reslo(k,*) = reslo(sel[k],*)/numslo[k]
endfor

remd = remd(0:nut-1,*)
reslo = reslo(0:nut-1,*)

  ; write data to the file
  name = '04'+string(mm, format='(I2)')+string(dd, format='(I2)')+'mdseeing2.dat'
openw, 1, name
print, 'Writing data to the file ',name
for k=0,nut-1 do begin
    printf, 1, rebinut[k], numslo[k], nummd[k],reslo(k,*), remd(k,*),$
    format='(F5.2,2I3,14F6.2)'
endfor
close, 1

!P.position = [0,0,0,0]
!P.charsize=1.8
!x.style=1.
!y.style=1.
!P.multi=[0,1,3]


!x.margin=[8,0.4]
!y.margin=[0,0]
!x.ticks=0
!x.tickinterval = 1
!x.minor = 1

xx = [utstart, utend]
yy = [0,1]
ylab = 0.8


k = 2
    plot, rebinut,  (remd(*,k)/6.8)^0.6  ,yr=yy
    oplot, rebinut, (reslo(*,k)/6.8)^0.6, psym=2
    xyouts, rebinut[2],ylab, 'Above 0.5km'    
    xyouts, 5,ylab,  string(dd, format='(I2)')+'/'+string(mm, format='(I2)')


k = 1
    plot, rebinut,  (remd(*,k)/6.8)^0.6  ,yr=yy
    oplot, rebinut, (reslo(*,k)/6.8)^0.6, psym=2
    xyouts, rebinut[2],ylab, 'At 0.5km'    

k = 0
    plot, rebinut,  (remd(*,k)/6.8)^0.6  ,yr=yy
    oplot, rebinut, (reslo(*,k)/6.8)^0.6, psym=2
    xyouts, rebinut[2],ylab, 'Below 0.5km'    

;stop
  


end