PRO cross_cats, ra1, dec1, ra2, dec2, radius, index, multi = multi, r_min = r_min
;
; PURPOSE:
; Identify overlap sources in two catalogs by searching source
; pairs with separations withing a certain value.
;
; INPUTS:
; ra1, dec1 - double precision vectors that contain source ra
; and dec in the first catalog, in unit of degree.
; ra2, dec2 - ra, dec of the second catalog
; radius - searching radius, in arcsec
; multi - set to 1 to avoid two sources having the same counterpart.
; In this case, the counterpart will be assigned to the source
; that is closest to the counterpart. Could be time consuming.
;
; OUTPUT:
; index - integer vector that has the same number of elements
; as the first catalog indicating which source in cat2
; is the counterpart. If -1, the corresponding cat1
; source doesn't have a counterpart in cat2.
;
; NOTE: 1. It doesn't work if your sources are close (<1 arcmin, maybe)
; to the poles.
; 2. The subroutine cross_dogs.pro is required.
;
;
If n_elements(multi) EQ 0 THEN multi = 0
n1 = n_elements(ra1)
IF n1 LE 2000 THEN cross_dogs, ra1, dec1, ra2, dec2, radius, index, multi = multi
IF n1 LE 2000 THEN goto, finish
index = lon64arr(n1)*0-1
; roughly 1000-2000 sources per block, if uniformly distributed
ncut = ceil(sqrt(ceil(n1/2000.0)))
ra_max = max(ra1)
ra_min = min(ra1)
dec_max = max(dec1)
dec_min = min(dec1)
dra = (ra_max - ra_min)/ncut
ddec = (dec_max - dec_min)/ncut
;order1 = sort(ra1)
;order2 = sort(dec1)
;dn = n1*1.0/ncut
FOR i=0,ncut-1 DO BEGIN
FOR j=0,ncut-1 DO BEGIN
A = where(ra1 GE i*dra+ra_min AND ra1 LE (i+1)*dra+ra_min AND $
dec1 GE j*ddec+dec_min AND dec1 LE (j+1)*ddec+dec_min)
IF total(A) NE -1 THEN BEGIN
cross_dogs, ra1[A], dec1[A], ra2, dec2, radius, index2, multi = multi
index[A] = index2
ENDIF
; why plan B doesn't work well??
; B = where(order1 GE i*dn AND order1 LE (i+1)*dn AND $
; order2 GE j*dn AND order2 LE (j+1)*dn)
; IF total(B) NE -1 THEN BEGIN
; cross_dogs, ra1[order1[B]], dec1[order2[B]], ra2, dec2, radius, index2, multi = multi
; index[B] = index2
; ENDIF
ENDFOR
ENDFOR
finish:
END