gtiseg.pro
4.8 KB
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
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
;+
; NAME:
; GTISEG
;
; AUTHOR:
; Craig B. Markwardt, NASA/GSFC Code 662, Greenbelt, MD 20770
; craigm@lheamail.gsfc.nasa.gov
;
; PURPOSE:
; Convert a list of times to a set of Good Time Intervals (GTIs)
;
; CALLING SEQUENCE:
; GTI = GTISEG(TIMES, COUNT=COUNT, INDICES=INDICES, $
; MAXGAP=MAXGAP, MINGTI=MINGTI)
;
; DESCRIPTION:
;
; The function GTISEG accepts an array of times and converts
; adjacent data into good time intervals (GTIs).
;
; Elements of the array are clustered into intervals based on the
; gaps between times. If the gaps are small enough then the times
; are grouped into a single interval. If a gap exceeds MAXGAP, then
; an interruption occurs and at least two intervals are formed.
; Thus, the keyword parameter MAXGAP essentially determines how many
; and where the intervals will be formed.
;
; If the time samples are regularly spaced -- aside from gaps --
; then MAXGAP should be set to a number slightly larger than the
; spacing to prevent roundoff errors. By default MAXGAP is set to
; the difference between the first and second samples.
;
; For GTISEG, the samples do not need to be regularly spaced, but
; they *must* be given in ascending order. Arrays can be sorted
; with the SORT function. The primary difference between GTISEG and
; MASK2GTI is that MASK2GTI assumes the time samples are regularly
; spaced while GTISEG does not. Also, MASK2GTI allows intervals to
; be enlarged or shrunk.
;
; It should be noted that this function is not constrained to
; operation only on time arrays. It should work on any
; one-dimensional quantity with intervals.
;
; INPUTS:
;
; TIME - an array of times in ascending order.
;
;
; KEYWORDS:
;
; MAXGAP - a scalar, the maximum gap between time samples before a
; new interval is created. Samples with gaps smaller than
; this value are grouped into a single GTI.
; Default: TIME(1) - TIME(0)
;
; MINGTI - the smallest possible GTI. Any interval smaller than
; MINGTI is discarded.
; Default: 0 (all intervals are accepted)
;
; COUNT - upon return, the number of resulting intervals. A value
; of zero indicates no good time intervals.
;
; INDICES - upon return, a 2xCOUNT array of integers which give the
; indices of samples which lie within each interval. The
; times TIME(INDICES(0,i) : INDICES(1,i)) fall within the
; ith interval.
;
; RETURNS:
;
; A new GTI array containing the enlarged or shrunken intervals.
; The array is 2xCOUNT where COUNT is the number of resulting
; intervals. GTI(*,i) represents the start and stop times of
; interval number i. The intervals are non-overlapping and
; time-ordered.
;
; If COUNT is zero then the returned array is a scalar value of
; zero, indicating no good intervals were found.
;
; SEE ALSO:
;
; MASK2GTI, GTITRIM, GTIMERGE, GTIWHERE
;
; MODIFICATION HISTORY:
; Written, CM, 1999-2001
; Documented, CM, Apr 2001
; MINGTI now works as documented, in that segments *equal* to MINGTI
; are now accepted, CM, 30 Oct 2007
; MINGTI now also affects INDICES, CM, 03 Mar 2008
; Handle case of scalar input, CM, 2016-04-15
;
; $Id: gtiseg.pro,v 1.8 2016/05/19 16:12:02 cmarkwar Exp $
;
;-
; Copyright (C) 1999-2001, 2007, 2008, 2016, Craig Markwardt
; This software is provided as is without any warranty whatsoever.
; Permission to use, copy, modify, and distribute modified or
; unmodified copies is granted, provided this copyright and disclaimer
; are included unchanged.
;-
function gtiseg, time, maxgap=maxgap, mingti=mingti, $
indices=indices, count=count, query=query
if keyword_set(query) then return, 1
count = 0L
if n_params() EQ 0 then begin
message, 'GTI = GTISEG(TIME, [COUNT=,] [INDICES=,] [MAXGAP=,] [MINGTI=])',$
/info
return, 0
endif
if n_elements(maxgap) EQ 0 then maxgap = time(1) - time(0)
nt = n_elements(time)
tdiff = [0d] & ntseg = 0L
if nt GT 1 then begin
tdiff = time(1:*) - time
whdiff = [-1L, where(tdiff GT maxgap(0), ntseg), nt-1]
endif
if (nt LE 1) OR (ntseg EQ 0) then whdiff = [-1L, nt-1]
ntseg = ntseg + 1
mintdiff = min(tdiff) > 0
indices = reform(lonarr(2, ntseg), 2, ntseg, /overwrite)
indices(0,*) = whdiff(0:ntseg-1)+1
indices(1,*) = whdiff(1:ntseg)
tgti = reform(make_array(2, ntseg, value=time(0)*0), 2, ntseg, /overwrite)
tgti(0,*) = time(indices(0,*))
tgti(1,*) = time(indices(1,*)) + mintdiff
if n_elements(mingti) GT 0 then begin
wh = where(tgti(1,*)-tgti(0,*) GE mingti(0), ntseg)
if ntseg GT 0 then begin
tgti = reform(tgti(*,wh),2,ntseg)
indices = reform(indices(*,wh),2,ntseg)
endif else begin
tgti = -1L
indices = -1L
endelse
endif
count = ntseg
return, tgti
end