pro vpd1,vx1,vy1,step,thr,im
;
; modification of vpd
; white background, dark coordinates
; step ... the arrow will be ploted for each step-th point
; thr  ... threshold to cut out too high values in vx1,vy1
; im   ... underlying image (e.g. avg. image of the movie)
;
vx=vx1>(-thr)
vx=vx<thr
vy=vy1>(-thr)
vy=vy<thr
;
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      ; differs from vpd

;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)
s=size(mag)
lmax=max(mag)
print,'Max and mean of velocity modulus:',lmax,mean(mag)
vect=fltarr(s(1),s(2))
mul=1.5*step/lmax		; 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=1,ystyle=1,$
/noerase,color=0,charsize=1 ; no color and charsize specification in vpd

for i=k,s(1)-k,k do begin
  for j=k,s(2)-k,k do begin
    x0=i
    x1=i+mul*vx(i,j)*2
    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))
      x3=x1-sin(a1-angle)*lte*vx(i,j)/abs(vx(i,j))
    y0=j
    y1=j+mul*vy(i,j)*2
      y2=y1-cos(a1+angle)*lte*vy(i,j)/abs(vy(i,j))
      y3=y1-cos(a1-angle)*lte*vy(i,j)/abs(vy(i,j))
;      plots,[x0,x1,x2,x1,x3],[y0,y1,y2,y1,y3],color=235,thick=4
     plots,[x0,x1,x2,x3,x1],[y0,y1,y2,y3,y1],color=0,thick=1
    endfor
endfor
;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