function two_src_visibility,vis_base,flux,area,x,y,xyoffset,PLOT=plot ; PURPOSE: returns complex visibility of two sources on RHESSI uvcircles as ; defined by vis_bas ; INPUTS: ; vis_base = RHESSI visibility structure (used only for [u,v] coordinates) ; Flux = float(2) = flux vector for two sources ; Area = float(2) = area vector for two sources ; = !pi*[a1^2,a2^2] where a1, a2 are Gaussian widths in asec ; x = float(2) = x coordinate of two sources relative to map center (asec) ; y = float(2) = y coordinate of two sources relative to map center (asec) ;OPTIONAL INPUTS: ; pa = vector of position angles (defaults to 32 values on 0 to 2*!pi) ; /plot to plot the amplitude profile ; ; OUTPUT: ; complex visibility of two sources ;;OPTIONAL INPUTS: ; ejs dec 12, 2006 ;pitch=(hsi_grid_parameters()).pitch ;kx=cos(pa) # (1/pitch) ;ky=sin(pa) # (1/pitch) u=vis_base.u v=vis_base.v pa=atan(v,u) k=sqrt(u^2+v^2) x1=x-xyoffset[0] ; shift to center of map y1=y-xyoffset[1] relvis1= exp(-2*!pi*k^2*area[0]) relvis2= exp(-2*!pi*k^2*area[1]) phz1 = 2*!pi*(u*x1[0]+v*y1[0]) phz2 = 2*!pi*(u*x1[1]+v*y1[1]) i=complex(0,1) ;stop vis_out=vis_base vis_out.obsvis = relvis1*flux[0]*exp(i*phz1) + relvis2*flux[1]*exp(i*phz2) amp=abs(Vis_out.obsvis) if keyword_set(plot) then begin window,1 ; plot,amp,/xst,/yst,tit='TWO-SOURCE VISIBILITY AMPLITUDE',/ylog,yra=[0.01,1]*total(flux),psym=-6,symsiz=0.3 plot,amp,/xst,/yst,tit='TWO-SOURCE VISIBILITY AMPLITUDE',yra=[0.01,1]*total(flux),psym=-6,symsiz=0.3 for j=0,8 do oplot,32*[j,j],[0.01,1]*total(flux) legend,/right,/bottom,['Flux='+strcompress(flux[0],/rem)+strcompress(flux[1]),$ 'xpos='+strcompress(x[0],/rem)+strcompress(x[1]),$ 'ypos='+strcompress(y[0],/rem)+strcompress(y[1]),$ 'area='+strcompress(area[0],/rem)+strcompress(area[1])],/clear,textcolor=1 endif ;stop return,vis_out end iplot=1 restore,'MEMfit_sigvis_phzradset_2003-06-17_225310-225440E30-60ff.sav',/ver ; get position angles: pa=(360+fix(atan(vis.v,vis.u)/!dtor)) mod 360 ; force PA to be in [0,360] scpa=vis.isc+pa/360. ; NOTE THIS IS DELIBERATLY NOT THE GH FORM WHICH GOES FROM 0 TO 180 amp=abs(vis.obsvis) s=sort(scpa) & sscpa=scpa(s) & samp=amp(s) sisc=(vis.isc)[s] ss=sisc(uniq(sisc)) & ws=[ss,n_elements(samp)-1] for j=min(ss),max(ss) do ws(j-1)=min(where(scpa ge j)) xt='SC+PA/360' if iplot ne 0 then begin window,0 plot,sscpa,samp,psym=-6,symsiz=0.3,tit='30617 visibility amplitudes',/xst,xtit=xt for j=0,8 do oplot,j*[1,1],[0,1000],line=2 legend,['date='+strcompress(anytim(vis[0].trange[0],/date,/ecs)),$ 'time='+strcompress(anytim(vis[0].trange[0],/ecs,/time)),$ 'energy band='+strcompress(fix(vis[0].erange[0]),/rem)+'-'+strcompress(fix(vis[0].erange[1]),/rem),$ 'phzrad=7'] write_png,'30617_vis_amplitudes_.png',255b-tvrd() endif fwhm=[4.7,3.8] ; from hsi_vis_fwdfit (chisq=1.55) a=fwhm/2.35 area=!pi*[a,a]^2 flux=[89.,74] x=[-797.89,-783.16] + 790.50 y=[-148.12,-128.93] + 138.00 vism=two_src_visibility(vis,flux,area,x,y,vis[0].xyoffset) pa=(360+fix(atan(vis.v,vis.u)/!dtor)) mod 360 ; force PA to be in [0,360] scpam=vism.isc+pa/360. ; NOTE THIS IS DELIBERATLY NOT THE GH FORM WHICH GOES FROM 0 TO 180 mamp=abs(vism.obsvis) s=sort(scpam) & sscpa=scpam(s) & smamp=mamp(s) siscm=(vism.isc)[s] ss=siscm(uniq(sisc)) & ws=[ss,n_elements(smamp)-1] for j=min(ss),max(ss) do ws(j-1)=min(where(scpam ge j)) xt='SC+PA/360' window,2 s1=findgen(n_elements(vis)) plot,sscpa,smamp,psym=-6,symsiz=0.9,tit='30617 VIS MODEL AMPLITUDES',/xst,xtit=xt legend,['fwhm='+strtrim(fwhm(0))+strtrim(fwhm[1]),$ 'xpos='+strtrim(x[0])+strtrim(x[1]),$ 'ypos='+strtrim(y[0])+strtrim(y[1]),$ 'flux='+strtrim(flux[0])+strtrim(flux[1])] wm=lonarr(9) for j=0,8 do wm(j)=min(where(vism.isc eq j)) for j=0,8 do oplot,j*[1,1],[0,1000],line=2 loadct,23 oplot,sscpa,samp,psym=-6,symsiz=0.3,color=250 legend,/bottom,/right,['2003/06/17 Amplitudes'],color=[254],box=0,position=[130,18] w=lonarr(9) for j=1,8 do oplot,w[j]*[1,1],[0,1000],line=1,color=250 loadct,0 legend,/bottom,/right,['Model Amplitudes'],box=0,position=[130,10] write_png,'30617_fwdfit_amplitudes_.png',255b-tvrd(/true) stop window,3 yra=[-150,150] plot,imaginary(vis.obsvis),psym=-6,symsiz=0.3,tit='30617 real & imaginary components',yra=yra,/xst oplot,float(vis.obsvis),line=1 legend,/right,/bott,['fwhm='+strtrim(fwhm(0))+strtrim(fwhm[1]),$ 'xpos='+strtrim(x[0])+strtrim(x[1]),$ 'ypos='+strtrim(y[0])+strtrim(y[1]),$ 'flux='+strtrim(flux[0])+strtrim(flux[1])] legend,['REAL(vis)','IMAG(vis)'],line=[1,0] for j=1,8 do oplot,w[j]*[1,1],[-1000,1000],line=2 write_png,'30617_real_imag_vis.png',255b-tvrd() window,4 plot,imaginary(vism.obsvis),psym=-6,symsiz=0.3,tit='VIS_FWDFIT MODEL real imaginary components',yra=yra,/xst oplot,float(vism.obsvis),line=1 for j=1,8 do oplot,wm[j]*[1,1],[-1000,1000],line=2 legend,['REAL(vis)','IMAG(vis)'],line=[1,0] write_png,'30617_vis_fwdfit_real_imag_vis.png',255b-tvrd() window,5 phase=atan(imaginary(vism.obsvis),float(vism.obsvis))/!dtor yra=[-180,180] plot,phase,psym=-6,symsiz=0.3,tit='MODEL PHASE (DEGREES)',/yst,/xst for j=1,8 do oplot,w[j]*[1,1],[-1000,1000],line=2 write_png,'model_phase.png',255b-tvrd() ;stop message,'Done' ;write_png,'2src_visF200-100W3-3.png',255b-tvrd() ;area=!pi*[3,3]^2 ;flux=[100.,100.] ;x=[-20,20] ;y=[-20,20] ;vis=two_src_visibility(flux,area,x,y,/plot) ;write_png,'2src_visF100-100W3-3.png',255b-tvrd() ;area=!pi*[2,3]^2 ;flux=[200.,100.] ;x=[-20,20] ;y=[-20,20] ;vis=two_src_visibility(flux,area,x,y,/plot) ;write_png,'2src_visF200-100W2-3.png',255b-tvrd() ;area=!pi*[2,4]^2 ;flux=[200.,100.] ;x=[-20,20] ;y=[-20,20] ;vis=two_src_visibility(flux,area,x,y,/plot) ;write_png,'2src_visF200-100W2-4.png',255b-tvrd() end