;+ ; Contains the dfe_main_full procedure ; ; :Author: ; Baptiste Cecconi ; ; :History: ; 2008/01/23: Created ; ; 2008/01/23: Last Edit ;- ; ;+ ; dfe_main_full is a procedure that <behavior desc here> ; ; :Uses: ; crossp1, read_antenna_set ; ; :Params: ; df_result: in, required, type=sometype ; A parameter named df_result ; data: in, required, type=sometype ; A parameter named data ; ephem: in, required, type=sometype ; A parameter named ephem ; antenna_file: in, required, type=sometype ; A parameter named antenna_file ; sn: in, required, type=sometype ; A parameter named sn ; nan_exception: in, required, type=sometype ; A parameter named nan_exception ; ; :Keywords: ; VERBOSE: in, optional, type=sometype ; A keyword named VERBOSE ; TEST: in, optional, type=sometype ; A keyword named TEST ;- PRO dfe_main_full,df_result,data,ephem,antenna_file,sn,nan_exception, $ VERBOSE=VERBOSE,TEST=TEST ; ----------------------------------------------------------------------------------- ; all angles are in radian ; ----------------------------------------------------------------------------------- if antenna_file eq 'calDec04' then begin read_antenna_set,antenna_HF,'calDec04_H12',/rad read_antenna_set,antenna_LF,'calDec04_ABC',/rad endif else begin read_antenna_set,antenna_HF,antenna_file,/rad antenna_LF = antenna_HF endelse f_antenna_lim = 320. ndim = n_elements(df_result) nan_exception = bytarr(ndim) ; -------------------------------------------------------------------- ; antenna parameters initialization ; -------------------------------------------------------------------- alZtmp = [antenna_LF.Z.al,antenna_HF.Z.al] alZ = alZtmp(data.f ge f_antenna_lim) beZtmp = [antenna_LF.Z.be,antenna_HF.Z.be] beZ = beZtmp(data.f ge f_antenna_lim) hZtmp = [antenna_LF.Z.h,antenna_HF.Z.h] hZ = hZtmp(data.f ge f_antenna_lim) alXtmp = [[0.,antenna_LF.Xp.al,antenna_LF.Xm.al,antenna_LF.Dip.al],[0.,antenna_HF.Xp.al,antenna_HF.Xm.al,antenna_HF.Dip.al]] alX = alXtmp(data.ant mod 10,data.f ge f_antenna_lim) beXtmp = [[0.,antenna_LF.Xp.be,antenna_LF.Xm.be,antenna_LF.Dip.be],[0.,antenna_HF.Xp.be,antenna_HF.Xm.be,antenna_HF.Dip.be]] beX = beXtmp(data.ant mod 10,data.f ge f_antenna_lim) hXtmp = [[0.,antenna_LF.Xp.h,antenna_LF.Xm.h,antenna_LF.Dip.h],[0.,antenna_HF.Xp.h,antenna_HF.Xm.h,antenna_HF.Dip.h]] hX = hXtmp(data.ant mod 10,data.f ge f_antenna_lim) ; --------------------------------------------------------------------------- ; Computing ephemeris in 2-antenna system frame ; NB: ; 2-antenna system frame = hX and hZ in xOy plane, y bissecting antenna directions ; --------------------------------------------------------------------------- vhZ = make_vect_sph(1.,alZ,beZ) vhX = make_vect_sph(1.,alX,beX) qhZ = q_vmake(vhZ) qhX = q_vmake(vhX) qsrc = q_vmake(make_vect_sph(1,ephem.th,ephem.ph)) e3 = [0.,0.,1.] vz1 = crossp1(vhz,vhx,/norm) vr1 = crossp1(vz1,e3,/norm) be1 = angular_distance(vz1,e3) qr1 = q_make(be1,vr1) qZ1 = q_rot(qr1,qhZ) qX1 = q_rot(qr1,qhX) phx1 = atan(qx1(2,*),qx1(1,*)) phZ1 = atan(qz1(2,*),qz1(1,*)) be2 = reform((!pi - (phZ1-phx1))/2.-phX1) qr2 = q_make(be2,e3) qz2 = q_rot(qr2,qz1) qx2 = q_rot(qr2,qx1) alX2 = acos(qX2(3,*)) alZ2 = acos(qZ2(3,*)) beX2 = atan(qX2(2,*),qX2(1,*)) beZ2 = atan(qZ2(2,*),qZ2(1,*)) qsrc2 = q_rot(qr2,q_rot(qr1,qsrc)) theph2 = acos(qsrc2(3,*)) pheph2 = atan(qsrc2(2,*),qsrc2(1,*)) ;!p.multi=[0,2,2] ;plot,alz2,psym=1,ytit='alz' ;plot,alx2,psym=1,ytit='alx' ;plot,bez2,psym=1,ytit='bez' ;plot,bex2,psym=1,ytit='bex' ;stop ; --------------------------------------------------------------------------- ; Preparing Useful Quantities... ; --------------------------------------------------------------------------- P = fltarr(ndim) & Q=P & R=P ws = where(data.autoZ ne 0. or data.autoX ne 0.) P(ws) = 2.*hX(ws)*hZ(ws)*data(ws).crossR/ $ (hx(ws)^2.*data(ws).autoZ+hZ(ws)^2.*data(ws).autoX) Q(ws) = (hx(ws)^2.*data(ws).autoZ-hZ(ws)^2.*data(ws).autoX)/ $ (hx(ws)^2.*data(ws).autoZ+hZ(ws)^2.*data(ws).autoX) R(ws) = 2.*hX(ws)*hZ(ws)*data(ws).crossI/ $ (hx(ws)^2.*data(ws).autoZ+hZ(ws)^2.*data(ws).autoX) c = sin(alX2)*sin(alZ2)*cos(beX2-beZ2)+cos(alX2)*cos(alZ2) Pc = (P-c)/(1.-P*c) Qc = Q*sqrt(1.-c^2.)/(1.-P*c) Rc = R*sqrt(1.-c^2.)/(1.-P*c) ; --------------------------------------------------------------------------- ; Determination of azimuth Ph ; --------------------------------------------------------------------------- ph2 = fltarr(ndim) ph2(ws) = 0.5*atan(-Qc(ws),Pc(ws)) ph2a = ph2 ph2b = ph2 ph2c = ph2+!pi ph2d = ph2+!pi ;test : ; IDL> plot,abs(Qc/Pc)<100,abs(tan(2.*pheph2))<100,psym=3 ; --------------------------------------------------------------------------- ; Determination of colatitude Th ; --------------------------------------------------------------------------- th2 = acos(sqrt((2./(1.+sqrt(Pc^2.+Qc^2.))-1.)>0.)<1.>(-1.)) th2a = th2 th2b = !pi-th2 th2c = th2 th2d = !pi-th2 ;test : ; IDL> plot,(2/(1+sqrt(Pc^2.+Qc^2.))-1),(cos(theph2)^2.),psym=3 ; --------------------------------------------------------------------------- ; Determination of flux S ; --------------------------------------------------------------------------- S = (data.autoX/hX^2.+data.autoZ/hZ^2.)*(1.+sqrt(Pc^2.+Qc^2.))/(1.+Pc*c) ; --------------------------------------------------------------------------- ; Determination of Circular polarization degree V ; --------------------------------------------------------------------------- V_a = Rc/(1.+sqrt(Pc^2.+Qc^2.))/cos(th2a) V_b = Rc/(1.+sqrt(Pc^2.+Qc^2.))/cos(th2b) V_c = Rc/(1.+sqrt(Pc^2.+Qc^2.))/cos(th2c) V_d = Rc/(1.+sqrt(Pc^2.+Qc^2.))/cos(th2d) ; --------------------------------------------------------------------------- ; Back into S/C frame from angles ; --------------------------------------------------------------------------- src2 = make_vect_sph(1.,th2a,ph2a) qsrc2 = q_vmake(src2) qsrc = q_rot(q_conj(qr1),q_rot(q_conj(qr2),qsrc2)) th_a = reform(acos(qsrc(3,*)<1.>(-1.))) ph_a = reform(atan(qsrc(2,*),qsrc(1,*))) mod (2.*!pi) src2 = make_vect_sph(1.,th2b,ph2b) qsrc2 = q_vmake(src2) qsrc = q_rot(q_conj(qr1),q_rot(q_conj(qr2),qsrc2)) th_b = reform(acos(qsrc(3,*)<1.>(-1.))) ph_b = reform(atan(qsrc(2,*),qsrc(1,*))) mod (2.*!pi) src2 = make_vect_sph(1.,th2c,ph2c) qsrc2 = q_vmake(src2) qsrc = q_rot(q_conj(qr1),q_rot(q_conj(qr2),qsrc2)) th_c = reform(acos(qsrc(3,*)<1.>(-1.))) ph_c = reform(atan(qsrc(2,*),qsrc(1,*))) mod (2.*!pi) src2 = make_vect_sph(1.,th2d,ph2d) qsrc2 = q_vmake(src2) qsrc = q_rot(q_conj(qr1),q_rot(q_conj(qr2),qsrc2)) th_d = reform(acos(qsrc(3,*)<1.>(-1.))) ph_d = reform(atan(qsrc(2,*),qsrc(1,*))) mod (2.*!pi) ; --------------------------------------------------------------------------- ; Filling data_n3e structure with DFe results ; --------------------------------------------------------------------------- v = transpose([[v_a],[v_b],[v_c],[v_d]]) th = transpose([[th_a],[th_b],[th_c],[th_d]]) ph = transpose([[ph_a],[ph_b],[ph_c],[ph_d]]) df_result.ydh = data.ydh df_result.num = data.num df_result.s = reform(s,ndim) df_result.v = reform(v,4,ndim) df_result.th = reform(th,4,ndim) df_result.ph = reform(ph,4,ndim) df_result.sn = sn return end