gtitrim.pro
6.14 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
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
;+
; NAME:
; GTITRIM
;
; AUTHOR:
; Craig B. Markwardt, NASA/GSFC Code 662, Greenbelt, MD 20770
; craigm@lheamail.gsfc.nasa.gov
;
; PURPOSE:
; Normalize a Good Time Interval (GTI) - no overlapping and adjoining
;
; CALLING SEQUENCE:
; NEWGTI = GTITRIM(GTI, COUNT=, MAXGAP=, MINGTI=)
;
; DESCRIPTION:
;
; A good time interval is by definition a list of intervals which
; represent "good" or acceptable portions of the real number line.
; In this library a GTI is a 2xNINTERVAL array where NINTERVAL is
; the number of intervals.
;
; The numbers in the array represent the start and stop times of
; each interval. Thus, the array [[0,10],[20,30],[40,50]] represent
; intervals ranging from 0-10, 20-30 and 40-50. Formally, this
; example GTI represents times which would satisfy the following
; expression, for each time T and interval number i:
; T GE GTI(0,i) AND T LT GTI(1,i)
; Note that the endpoint is closed on the left but open on the
; right.
;
; However, not every 2xNINTERVAL array is a valid or "normalized"
; GTI as used by this library. The array must satisfy several
; conditions:
; * time ordered (ascending)
; * no overlapping intervals
; * no adjoining intervals (intervals that start and stop at the
; same point; e.g. the point 10 in this array [[0,10],[10,20]])
;
; A user who desires to create his or her own GTI array can proceed
; as follows.
;
; First, the array is placed in time order. This can be
; accomplished simply using the built-in function SORT. This
; statement sorts the array by start times.
;
; GTI = GTI(*, SORT(GTI(0,*)))
;
; Second, the GTITRIM function is used to fully normalize the set of
; intervals:
;
; GTI = GTITRIM(GTI)
;
; After these two procedures the GTI is considered valid and can be
; passed to the other routines of the library. Of course if the
; user can guarantee the above requirements without using GTITRIM
; then this is acceptable as well.
;
; 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:
;
; GTI - a 2xNINTERVAL array where NINTERVAL is the number of
; intervals. GTI(*,i) represents the start and stop times of
; interval number i. The intervals must be non-overlapping
; and time-ordered (use GTITRIM to achieve this).
;
; A scalar value of zero indicates that the GTI is empty, ie,
; there are no good intervals.
;
; KEYWORDS:
;
; MAXGAP - Maximum allowable gap for merging existing good time
; intervals. Intervals with gaps smaller than MAXGAP will
; be combined into a single interval.
; Default: 0 (any gap keeps intervals separate)
;
; MINGTI - Minimum size interval. If any interval is smaller than
; MINGTI then it is discarded.
; Default: 0 (all intervals are preserved)
;
; COUNT - upon return, the number of resulting intervals. A value
; of zero indicates no good time intervals.
;
;
; RETURNS:
;
; A new GTI array containing the normalized 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:
;
; GTIMERGE
;
; MODIFICATION HISTORY:
; Written, CM, 1997-2001
; Documented, CM, Apr 2001
; Corrected bug which bypassed MIN/MAXGTI, CM, 20 Jul 2003
;
; $Id: gtitrim.pro,v 1.7 2006/10/22 09:50:09 craigm Exp $
;
;-
; Copyright (C) 1997-2001, 2003, 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 gtitrim, gti, maxgap=maxgap, mingti=mingti, maxgti=maxgti, $
count=count, query=query
if keyword_set(query) then return, 1
count = 0L
ngti = n_elements(gti)
if n_elements(maxgap) EQ 0 then maxgap = 0D
if n_elements(mingti) EQ 0 then mingti = 0D
;; Case of empty GTI
if ngti LT 2 then return, 0L
;; Case of single GTI
if ngti EQ 2 then begin
if gti(1) - gti(0) LT mingti then return, 0L
count = 1L
newgti = reform(gti, 2, 1)
goto, GTI_CHECKS
endif
newgti = gti
ngti = ngti/2
;; Check for gaps, spaces between GTIs
gaps = newgti(0,1:*) - newgti(1,*)
wh = where(gaps GT maxgap, ct)
if ct EQ 0 then begin
;; All of the gaps are too small
newgti = reform([min(newgti), max(newgti)], 2, 1)
endif else if ct GT 0 AND ct LT ngti-1 then begin
;; Remake the GTI, removing the gaps
vgti = reform(make_array(2, ct+1, value=newgti(0)), 2, ct+1, /overwrite)
vgti(0,0) = min(newgti)
vgti(1,ct) = max(newgti)
vgti(1,0:ct-1) = newgti(1,wh)
vgti(0,1:ct) = newgti(0,wh+1)
newgti = temporary(vgti)
endif
GTI_CHECKS:
;; Remove too-small GTIs
dur = newgti(1,*) - newgti(0,*)
wh = where(dur GE mingti, count)
if count GT 0 then newgti = newgti(*,wh) else newgti = 0L
if count GT 0 then newgti = reform(newgti, 2, count, /overwrite)
;;
if count GT 0 AND n_elements(maxgti) GT 0 then begin
maxgti1 = maxgti(0)
nper = (newgti(1,*)-newgti(0,*)) / maxgti1
iper = ceil(nper)
totgti = total(iper)
if maxgti1 GT 0 AND totgti GT count then begin
vzero = newgti(0) & vzero(0) = 0
newgti2 = make_array(2, totgti, value=vzero)
j = 0L
for i = 0L, count-1 do begin
if iper(i) EQ 1 then begin
newgti2(*,j) = newgti(*,i)
endif else begin
for k = 0, iper(i)-1 do $
newgti2(*,j+k) = [ newgti(0,i)+maxgti1*k, $
newgti(1,i)<(newgti(0,i)+maxgti1*(k+1))]
endelse
j = j + iper(i)
endfor
count = totgti
newgti = reform(newgti2,2,totgti)
endif
endif
return, newgti
end