pro INIT_DF, fichier,t97,co,az,com,azm,os,pos = pos ,interval = interval, $
path_ephatt = path_ephatt, xyz = xyz, km = km, rad = rad, quiet=quiet
if keyword_set(quiet) then quiet=1b else quiet=0b
if n_elements(os) eq 0 then os = 0
if not keyword_set(path_ephatt) then path_ephatt = getenv('ROOT_RPWS') + '/ephatt/''
fichier=strtrim(string(fichier),2)
READ_EPH, fichier, x, y, z, r, t97, rp, os, path=path_ephatt
xyz = [[x],[y],[z]]
if keyword_set(km) then xyz = xyz*Rp
READ_ATT, fichier, m, t97a, os, path=path_ephatt
n=n_elements(t97)
mm=fltarr(3,3,n)
for i=0,2 do for j=0,2 do mm(i,j,*)=interpol(m(i,j,*),t97a,t97)
m=mm
if keyword_set(pos) then ref=pos else ref=[0.,0.,0.]
if keyword_set(km) then ref = ref/Rp
nref = n_elements(ref)/3
ref = reform(ref,3,nref)
co=fltarr(n,nref)-999. & az=co
for j=0,nref-1 do begin
for i=0,n-1 do begin
vect=[-x(i)+ref(0,j),-y(i)+ref(1,j),-z(i)+ref(2,j)]
rr=sqrt(vect(0)^2+vect(1)^2+vect(2)^2)
vect=vect/rr
mm=reform(m(*,*,i))
res=mm#vect
res=long(res*100000l)/100000.
if abs(res(2)) lt 1. then begin
colat=acos(res(2))
if abs(res(0)/sin(colat)) le 1. then begin
azi=acos(res(0)/sin(colat))*!radeg
if res(1) lt 0. then azi=360.-azi
az(i,j)=azi
endif
colat=colat*!radeg
co(i,j)=colat
endif else if abs(res(2)) eq 1. then begin
co(i,j) = 0.
az(i,j) = 0.
endif
endfor
endfor
for j=0,nref-1 do begin
ok=where(co(*,j) ne -999 and az(*,j) ne -999)
t97=t97(ok) & & co=co(ok,*) & az=az(ok,*) & xyz = xyz(ok,*)
endfor
if keyword_set(interval) then begin
print," restricting time interval..."
test = where((t97 mod 1)*24. le interval(1) $
and (t97 mod 1)*24. ge interval(0))
t97 = t97(test)
co = co(test,*)
az = az(test,*)
endif
com = fltarr(nref) & azm=com
for j = 0,nref-1 do begin
if not quiet then print,'source #',j
s=stdev(co(*,j),comtmp) & com(j)=comtmp
if not quiet then print,'< colatitude > = ',comtmp,' +/- ',s,' (degres - 1 sigma)'
s=stdev(az(*,j),azmtmp) & azm(j)=azmtmp
if not quiet then print,' < azimuth > = ',azmtmp,' +/- ',s,' (degres - 1 sigma)'
endfor
if nref eq 1 then begin
co = reform(co)
az = reform(az)
com = com(0)
azm = azm(0)
endif
help,t97,co,az
if keyword_set(rad) then begin
co = co/!radeg
az = az/!radeg
com = com/!radeg
azm = azm/!radeg
endif
return
end