PRO OTROK3,inout,result ;Evaluation of results of OTRACK revised by OTRIM (input arrays INOUT). ;The mean velocity is determined by linear least squares fit. ;Output: array RESULT (12,n) Each row corresponds to one object. ; 0 - lifetime in minutes ; 1 - average position x (pixel coords) ; 2 - average position y ; 3 - mean velocity in x (km/s) ; 4 - stdev of velocity in x (km/s) ; 5 - mean velocity in y (km/s) ; 6 - stdev of velocity in y (km/s) ; 7 - mean velocity modulus (km/s) ; 8 - average peak intensity (units of Iph) ; 9 - stdev of intensity ; 10 - average effective diameter (") ; 11 - stdev of effective diameter ; ; 15 July 1996, Michal. ;************** PARAMETERS **************************** dtim=20 ;time difference between frames (s) scal=0.062 ;scale arcsec/pix sx =180 ;x-size of frames processed by OTRACK ;****************************************************** kms=scal*725. ;conversion pixels to km s=size(inout) n=s(3) ;number of objects result=fltarr(12,n) ;Loop over objects for i=0,n-1 do begin w=where(inout(0,*,i),m) ;subscripts of existing data, ;m - number of frames where object exists ;Lifetime result(0,i)=(m-1)*dtim/60. ;in minutes ;Average position pos=inout(1,w,i) ;1-D coordinate x=pos mod sx ;2-D coords y=pos/sx result(1,i)=round(total(x)/m) result(2,i)=round(total(y)/m) ;Average velocity ;Linear least squares fit (for formulae see Rektorys 1973, p.1065) ;Mean velocity is defined as a slope of x(t) and y(t) t=dindgen(m)*dtim ;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/m ssy=suty-sut*suy/m sst=sut2-sut*sut/m vx=ssx/sst ;slopes (velocity pix/s) vy=ssy/sst erx=sqrt((sux2-sux*sux/m-vx*ssx)/(m-2)/sst) ;sigma errors ery=sqrt((suy2-suy*suy/m-vy*ssy)/(m-2)/sst) result(3,i)=float(vx*kms) ;velocity x result(4,i)=float(erx*kms) ;error velocity x result(5,i)=float(vy*kms) ;velocity y result(6,i)=float(ery*kms) ;error velocity y result(7,i)=sqrt(result(3,i)^2+result(5,i)^2) ;velocity modulus ;Peak intensity, mean, stdev result(9,i)=stdev(inout(0,w,i)/10000.,mea) result(8,i)=mea ;Effective diameter, mean, stdev diam=scal*sqrt(4*inout(2,w,i)/!pi) result(11,i)=stdev(diam,mea) result(10,i)=mea endfor print,n,' objects' end