forked from low-sky/cprops
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathderiv_decimate_kernel.pro
113 lines (97 loc) · 3.7 KB
/
deriv_decimate_kernel.pro
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
function deriv_decimate_kernel, cube, merger, lmaxin, levels = levels $
, sigdiscont = sigdiscont, nredun = nredun $
, fscale = fscale
;+
; NAME:
; DERIV_DECIMATE_KERNEL
; PURPOSE:
; To decimate a kernel set based on the jumps in the moments of the
; emission before and after merging.
;
; CALLING SEQUENCE:
; new_kernels = DERIV_DECIMATE_KERNEL( cube, merger_matrix, kernels
; [ levels = levels, sigdistcont = sigdistcont, $
; nredun = nredun, fscale = fscale, noextrap = noextrap])
;
; INPUTS:
; CUBE -- A data cube
; MERGER_MATRIX -- A matrix generated by the MERGEFIND program that
; tracks the levels where various kernels merge.
; KERNELS -- The indices of local maxima within CUBE
;
; KEYWORD PARAMETERS:
; LEVELS -- The levels with which to contour the data.
; SIGDISCONT -- The fractional change in the moments to register a
; discontinuity, a scalar or a 4-elt. vector
; corresponding to the jumps in [sigma_x, sigma_y,
; sigma_v and int_flux]
; NREDUN -- The number of moments that need to jump to count the
; kernels as being discontinuous.
; FSCALE -- The scale factor to multiply the fractional change in
; integrated flux by to compare with the same value for
; the second moments.
;
; OUTPUTS:
;
;
; MODIFICATION HISTORY:
;
; Tue Jun 18 21:05:40 2013, erosolo <erosolo@>
;
; Added ability to actually skip derivative decimation for
; performance in large data cubes.
;
; Documented -- Fri Sep 2 16:20:12 2005, Erik Rosolowsky
;-
if total(sigdiscont) eq 0 then begin
message,'Skipping Derivative Decimation',/con
return, lmaxin
endif
vectorify, cube, mask = (cube eq cube) $
, x = x, y = y, v = v, t = t, id = id $
, sz = sz, indvec = indvec
lmaxes = lonarr(n_elements(lmaxin))
for k = 0, n_elements(lmaxin)-1 do $
lmaxes[k] = where(lmaxin[k] eq indvec)
xgrid = (intarr(n_elements(lmaxes))+1) ## indgen(n_elements(lmaxes))
ygrid = indgen(n_elements(lmaxes)) ## (intarr(n_elements(lmaxes))+1)
mask = xgrid gt ygrid
uselmax = lonarr(n_elements(lmaxes)) + 1
ind_of_interest = where(mask eq 1, nindofint)
while (nindofint gt 0) do begin
highest_merge = max(merger*mask, /NAN, maxind)
kernel1_ind = xgrid[maxind]
kernel2_ind = ygrid[maxind]
;
kernels = [lmaxes[kernel1_ind], lmaxes[kernel2_ind]]
; CHECK IF THE TWO KERNELS ARE DISCONTINUOUS IN PROPERTIES AT THE
; MERGER. IF SO, WE WILL KEEP BOTH.
keep = are_kernels_disc(x, y, v, t, kernels $
, shared_contour = highest_merge $
, levels = levels $
, sigdiscont = sigdiscont, nredun = nredun $
, fscale = fscale)
; IF NOT, WE KEEP ONLY THE HIGHER OF THE TWO KERNELS
if (keep eq 0) then begin
if (merger[kernel1_ind, kernel1_ind] gt merger[kernel2_ind, kernel2_ind]) then begin
uselmax[kernel2_ind] = 0
mask[kernel2_ind, *] = 0
mask[*, kernel2_ind] = 0
endif else begin
uselmax[kernel1_ind] = 0
mask[kernel1_ind, *] = 0
mask[*, kernel1_ind] = 0
endelse
endif else begin
mask[kernel1_ind, kernel2_ind] = 0
mask[kernel2_ind, kernel1_ind] = 0
endelse
ind_of_interest = where(mask eq 1, nindofint)
endwhile
; RETURN THE SURVIVING LOCAL MAXIMA
newlmax = lmaxin[where(uselmax eq 1)]
kernel_ct = strcompress(string(n_elements(newlmax)))
message, /con, 'Number of Kernels after derivative decimation: '+kernel_ct
return, newlmax
end ; OF DERIV_DECIMATE_KERNELS