PRO mantrack,mov,res ;Manual tracking of one object. The object has to be marked ;by the cursor in each frame. ;Input - segmented movie MOV ;Output - 2-D array RES: one row per frame, ; column 0: max. intensity of the object ; column 1: x-coordinate of the maximum ; column 2: y-coordinate of the maximum ; column 3: avg. intensity of the object ; column 4: d_eff of the object (") ; column 5: instantaneous velocity (km/s) on_error,1 ;************** PARAMETERS **************************** lag=20 ;time difference between frames (s) scal=0.062 ;scale arcsec/pix ;****************************************************** si=size(mov) if si(0) ne 3 then message,'*** Only movies are allowed!' scale=scal*725. ;km/pix loadct,3 window,2,xsiz=si(1),ysiz=si(2),title='frame',xpos=100,ypos=600 window,3,xsiz=si(1),ysiz=si(2),title='object',xpos=120+si(1),ypos=600 res=fltarr(6,si(3)) for i=0,si(3)-1 do begin ima=reform(mov(*,*,i)) wima=bytarr(si(1),si(2)) OBJEKT1,ima,tab wset,2 tvscl,ima vedle: print,i,' Select object by cursor (click left)' cursor,x,y,/device,/down if !err eq 4 then goto,fin pos=y*si(1)+x w=where(tab(1,*) eq pos,zz) if zz eq 0 then goto,vedle objno=tab(0,w) ;find object obj=tab(1,where(tab(0,*) eq objno(0),c)) ;array of 1-D coords res(0,i)=max(ima(obj),whmax) ;intensity maximum maxpos=obj(whmax) ;1-D coord. of maximum res(1,i)=maxpos mod si(1) ;x-coord. of maximum res(2,i)=maxpos/si(1) ;y-coord. of maximum res(3,i)=total(ima(obj))/c ;average intensity res(4,i)=c ;area of object (pix) wima(obj)=1 wima(res(1,i)-1:res(1,i)+1,res(2,i))=0 wima(res(1,i),res(2,i)-1:res(2,i)+1)=0 wset,3 tvscl,wima endfor fin: res=res(*,where(res(0,*))) ;truncation sir=size(res) print,'Number of measured positions:',sir(2) ;normalization res(0,*)=res(0,*)/10000. ;intensity I/Iph res(3,*)=res(3,*)/10000. res(4,*)=scal*sqrt(4*res(4,*)/!pi) ;d_eff ;Average and instantaneous velocities window,0,xsiz=512,ysiz=382 window,1,xsiz=512 x=reform(res(1,*)) y=reform(res(2,*)) ;Average speed t=dindgen(sir(2))*lag ;time coordinate in seconds sut=total(t) ;linear least squares fit... sut2=total(t*t) sux=total(x) sux2=total(long(x)*x) suy=total(y) suy2=total(long(y)*y) sutx=total(t*x) suty=total(t*y) ssx=sutx-sut*sux/sir(2) ssy=suty-sut*suy/sir(2) sst=sut2-sut*sut/sir(2) vx=ssx/sst ;slopes (velocity pix/s) vy=ssy/sst erx=sqrt((sux2-sux*sux/sir(2)-vx*ssx)/(sir(2)-2)/sst) ery=sqrt((suy2-suy*suy/sir(2)-vy*ssy)/(sir(2)-2)/sst) erx=scale*erx ;sigma errors (km/s) ery=scale*ery v=scale*sqrt(vx*vx+vy*vy) ;velocity modulus print,'Average speed (km/s):',v print,'Errors (vx,vy):',erx,ery ;Interactive smoothing of trajectory and calculation of speed ;Initial parameters for SMSPLINE tens=5 ;tension (decreases with value) nods=3 ;number of nodes changepar: wset,1 !p.multi=[0,0,2,0,0] plot,x,ytitle='x-trajectory',/yst xs=SMSPLINE(x,tens,nods) oplot,xs,col=150 plot,y,ytitle='y-trajectory',/yst ys=SMSPLINE(y,tens,nods) oplot,ys,col=150 !p.multi=0 ;Speed vx=(xs-shift(xs,1))*scale/lag vy=(ys-shift(ys,1))*scale/lag vx=vx(1:sir(2)-1) vy=vy(1:sir(2)-1) v=sqrt(vx*vx+vy*vy) wset,0 plot,v,ytitle='velocity magnitude (km/s)' print,'Smoothing parameter =',nods read,'New smoothing parameter (ge 3) or accept (0): ',nods ; set new nods or finish if nods ge 3 then goto,changepar res(5,*)=[v(0),v] wdelete,1,2,3 end