pro vpd4,vx,vy,step,arlen,im,scale,WHITE=white ; ; version 4 with a possibility to scale the axes to arcseconds ; a bug was removed and it is not necessary to run it twice ; more reasonable scaling of arrow lengths, given by ARLEN ; modification of vpd ; white background, dark coordinates ; step ... the arrow will be ploted for each step-th point ; arlen... scaling of arrow length ; im ... underlying image (e.g. avg. image of the movie) ; scale... arcsec/pix scaling ; ; velocity unit = 2*step*arlen pixels ; if n_params() lt 6 then scale=1. plot,[0,1],/nodata,xrange=[0,1],yrange=[0,1],$ xstyle=5,ystyle=5,back=255,color=255,charsize=1.8,$ xtitle='arcsec',ytitle='arcsec' ;preparation, orig.charsize=1.5 px=!x.window*!d.x_vsize py=!y.window*!d.y_vsize sx=px(1)-px(0)+1 sy=py(1)-py(0)+1 print,'new image size:',sx,sy ;erase ;tv,bytarr(!d.x_vsize,!d.y_vsize)+255 ; white background ;set_plot,'ps' ;device,filename='d;\data\palma\map.ps',/color ;loadct,3 tvscl,congrid(bytscl(im),sx,sy,/cubic),px(0),py(0) ;tv,congrid(bytscl(im),sx,sy,/cubic),0,0,xsize=12.7,ysize=12.7,/centimeters ; mag=sqrt(vx^2+vy^2) ; unit - pixels/lag s=size(mag) print,'Max and mean of velocity modulus:',max(mag),mean(mag) vect=fltarr(s(1),s(2)) mul=arlen*step ; arrow length... k=fix(step) x0=0 x1=0 y0=0 y1=0 ;!p.position=[0,0,0.7155,1] ;plot,[x0,x1],[y0,y1],/nodata,xrange=[0,s(1)],yrange=[0,s(2)],title='vel',$ ; xstyle=1,ystyle=1,/noerase,/normal,xthick=5,ythick=5,charthick=5 plot,[x0,x1],/nodata,xrange=[0,s(1)],yrange=[0,s(2)],$ xstyle=5,ystyle=5,/noerase,color=0,charsize=1.8 for i=k,s(1)-k,k do begin for j=k,s(2)-k,k do begin x0=i ;beginning of arrow (x) x1=i+mul*vx(i,j)*2 ;end of arrow (x) angle=22.5*!dtor a1=abs(atan(vx(i,j)/vy(i,j))) st=mul*mag(i,j)*sin(angle) ct=mul*mag(i,j)*cos(angle) lte=sqrt(st^2+ct^2)/2. x2=x1-sin(angle+a1)*lte*vx(i,j)/abs(vx(i,j)) ;arrow head x3=x1-sin(a1-angle)*lte*vx(i,j)/abs(vx(i,j)) y0=j ;beginning of arrow (y) y1=j+mul*vy(i,j)*2 ;end of arrow (y) y2=y1-cos(a1+angle)*lte*vy(i,j)/abs(vy(i,j)) ;arrow head y3=y1-cos(a1-angle)*lte*vy(i,j)/abs(vy(i,j)) if keyword_set(white) then $ plots,[x0,x1,x2,x1,x3],[y0,y1,y2,y1,y3],color=255,thick=1 $ else $ plots,[x0,x1,x2,x1,x3],[y0,y1,y2,y1,y3],color=0,thick=1 ;open head ; plots,[x0,x1,x2,x3,x1],[y0,y1,y2,y3,y1],color=0,thick=1 ;closed head endfor endfor arrowunit=2*step*arlen plots,[0,arrowunit],[0,0],color=0,thick=4 ;velocity unit bar print,'arrow length unit =',arrowunit plot,[x0,x1],/nodata,xrange=[0,s(1)*scale],yrange=[0,s(2)*scale],$ xstyle=1,ystyle=1,/noerase,color=0,charsize=1.8,ticklen=(-0.02),$ xtitle='x [arcsec]',ytitle='y [arcsec]' ;box of axes ;div=divergen(vx,vy) ;lh=max(div)/2. ;l1=lh/5. ;contour,div,xstyle=5,ystyle=5,/noerase,$ ; nlevels=10,color=235,/normal,c_thick=2 ;device,/close ;set_plot,'x' end