PRO FLOWMAKR,MOV,LAG,FWHM,REB,VX,VY,BOXCAR=boxcar,ADIF=adif,CORR=corr, $
	QFIT2=qfit2,CROSSD=crossd
;+
; NAME:
;	FLOWMAKR (adapted FLOW_MAKER for array)
;
; PURPOSE:
;	Compute flow maps.
;
; INPUTS:
;	MOV = 3-D array containing the series of images
;	
;	LAG = time-lag between "reference" and "life" subseries
;	       (in number of frames)
;
;	FWHM = FWHM for smoothing window in pixels
;
;	REB = rebinning factor to change scale/range of Novembers's
;	method. First, MOV(*,*,i) is rebinned down by factor REB,
;	then flow-map is computed, and finally it is rebinned up
;	(interpolated) by factor REB to original size.
;
; KEYWORDS:
;	BOXCAR = if set, a boxcar window of width "FWHM" is used. Hence,
;		FWHM must be an odd number.
;
;	ADIF = uses an absolute differences algorithm.
;
;	CORR = uses a multiplicative algorithm. Default is the sum of
;		square of the local differences.
;
;	QFIT2 = uses 9 points fitting procedure.
;
;	CROSSD = uses cross derivative interpolation formulae.
;
; OUTPUTS:
;	VX,VY = X and Y components for the proper motion map.
;
; SIDE EFFECTS:
;	None.
;
; COMMON BLOCKS:
;	None.
;
; PROCEDURE:
;	It uses November's method of shifting BOTH images. The defaults for
;	the different methods are square differences for the matching, 
;	gaussian window for the smoothing and FIVEPOINT for the subpixel
;	extrem finding procedure.
;
; EXAMPLE:
;	IDL> FLOWMAKR,movie,2,8,2,vx,vy
;
;	Lag is 2 frames.
;	FWHM of gaussian window is 8 pixels.
;	Rebinning factor is 2.
;
; REFERENCES:
;	November,L.J. and Simon,G.W.: 1988, Ap.J., 333, 427
;	Darvann,T.: 1991, Master's Thesis, University of Oslo.
;
; MODIFICATION HISTORY:
;	Written by Roberto Luis Molowny Horas, May 1992.
;	Keywords added in July 1992.
;	FITS routines added in November 1992.
;	Modification for array input and rebinning
;	 - JAB, JKH & SOB, Oct 1994.
;
;-
ON_ERROR,2

; --------------- TUNABLE PARAMETERS -----------------------------

smo=7 ; 12  ; smoothing parameter for unsharp masking (pixels)
umq=0.6	      ; unsharp masking coefficient (0...no u.m., 1...full u.m.)
shf=1  ; shift in the correlation to evaluate the displacement map (pix)

; ----------------------------------------------------------------

	std1=fwhm*0.424661          ; conversion fwhm/std for Gaussian
	std=std1/reb

	s1 = SIZE(mov)		    ;Acquiring array dimensions.
	if s1(0) ne 3 then begin
		print,'Array must be 3-dimensional! Stopped.'
		stop
	endif

; -----------------------------------------------------------------

; Equalized image (unsharp masking and log scaling)...

mov1=mov

for i=0,s1(3)-1 do begin

  ima=SMOOTHE(mov(*,*,i),smo)
  ima=mov(*,*,i)-umq*ima

; Cut negative values (and make the logarithm) to enhance faint features

;  ------------ (other possibility - beg.)
;  ima=ima>0
;  (Scaling images to be edible for ALOG)
;  ima=100.*float(ima)/mean(ima)>1.  ; here is applied a selection that a BF
                                     ;  must be brighter than 1% of the mean
;  ------------ (other possibility - end)

; ima=ima>1		; careful with the minimum value - here
			;  we have the normalization at 10000 !
; ima=alog(ima)
  mov1(0,0,i)=ima

endfor

; ---------------------------------------------------------------

	s=[s1(0),fix(s1(1)/reb),fix(s1(2)/reb)]  ; rebinned size

	n=s1(3)-lag		; no of frames in subseries
	np=long(s(1))*long(s(2)) ; n_elements in a rebinned frame

	cc = FLTARR(s(1),s(2),3,3)	;The cum. correlation function.

	FOR k = 0,n-1 DO BEGIN
		a = congrid(mov1(*,*,k),s(1),s(2))
		b = congrid(mov1(*,*,k+lag),s(1),s(2))
		a = a - TOTAL(a)/np		;Remove mean.
		b = b - TOTAL(b)/np

		FOR i = -1,1 DO FOR j = -1,1 DO CASE 1 OF	;Methods.
			KEYWORD_SET(adif): BEGIN	;Absolute differences.
				cc(0,0,i+1,j+1) = cc(*,*,i+1,j+1)+$
				ABS(SHIFT(a,i*shf,j*shf)-SHIFT(b,-i*shf,-j*shf))
				END
			KEYWORD_SET(corr): BEGIN	;Cross products.
				cc(0,0,i+1,j+1) = cc(*,*,i+1,j+1)+$
				SHIFT(a,i*shf,j*shf) * SHIFT(b,-i*shf,-j*shf)
				END
			ELSE: BEGIN			;Square differences.
			dumb = SHIFT(a,i*shf,j*shf) - SHIFT(b,-i*shf,-j*shf)
				cc(0,0,i+1,j+1) = cc(*,*,i+1,j+1) + dumb*dumb
				dumb = 0	;This is faster than (...)^2
				END
			ENDCASE
		a = 0 & b = 0
	ENDFOR

	cc(0,0,0,0) = cc(1,*,*,*)		;Takes care of the edges.
	cc(0,0,0,0) = cc(*,1,*,*)
	cc(s(1)-1,0,0,0) = cc(s(1)-2,*,*,*)
	cc(0,s(2)-1,0,0) = cc(*,s(2)-2,*,*)

	FOR i = 0,2 DO FOR j = 0,2 DO IF KEYWORD_SET(boxcar) THEN $
	cc(0,0,i,j) = SMOOTHE(cc(*,*,i,j),fix(fwhm/reb)) ELSE $	;Boxcar...
	cc(0,0,i,j) = SCONVOL(cc(*,*,i,j),std=std)		;or gausian.

	CASE 1 OF
		KEYWORD_SET(qfit2): QFIT2,cc,vx,vy		;9-p. fitting,
		KEYWORD_SET(crossd): CROSSD,cc,vx,vy		;cross derivat.
		ELSE: FIVEPOIN,cc,vx,vy				;or default.
	ENDCASE
	cc = 0
	vx = 2.*shf*vx & vy = 2.*shf*vy		;Scales the result.

	vx=congrid(vx,s1(1),s1(2))*reb		;Restores original size
	vy=congrid(vy,s1(1),s1(2))*reb		; and values of flow map.
END