Chris@2
|
1 /*
|
Chris@2
|
2 Copyright (C) 2001, 2006 by Simon Dixon
|
Chris@2
|
3
|
Chris@2
|
4 This program is free software; you can redistribute it and/or modify
|
Chris@2
|
5 it under the terms of the GNU General Public License as published by
|
Chris@2
|
6 the Free Software Foundation; either version 2 of the License, or
|
Chris@2
|
7 (at your option) any later version.
|
Chris@2
|
8
|
Chris@2
|
9 This program is distributed in the hope that it will be useful,
|
Chris@2
|
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
|
Chris@2
|
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
Chris@2
|
12 GNU General Public License for more details.
|
Chris@2
|
13
|
Chris@2
|
14 You should have received a copy of the GNU General Public License along
|
Chris@2
|
15 with this program (the file gpl.txt); if not, download it from
|
Chris@2
|
16 http://www.gnu.org/licenses/gpl.txt or write to the
|
Chris@2
|
17 Free Software Foundation, Inc.,
|
Chris@2
|
18 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
|
Chris@2
|
19 */
|
Chris@2
|
20
|
Chris@2
|
21 package at.ofai.music.util;
|
Chris@2
|
22
|
Chris@2
|
23 import java.util.LinkedList;
|
Chris@2
|
24
|
Chris@2
|
25 public class Peaks {
|
Chris@2
|
26
|
Chris@2
|
27 public static boolean debug = false;
|
Chris@2
|
28 public static int pre = 3;
|
Chris@2
|
29 public static int post = 1;
|
Chris@2
|
30
|
Chris@2
|
31 /** General peak picking method for finding n local maxima in an array
|
Chris@2
|
32 * @param data input data
|
Chris@2
|
33 * @param peaks list of peak indexes
|
Chris@2
|
34 * @param width minimum distance between peaks
|
Chris@2
|
35 */
|
Chris@2
|
36 public static int findPeaks(double[] data, int[] peaks, int width) {
|
Chris@2
|
37 int peakCount = 0;
|
Chris@2
|
38 int maxp = 0;
|
Chris@2
|
39 int mid = 0;
|
Chris@2
|
40 int end = data.length;
|
Chris@2
|
41 while (mid < end) {
|
Chris@2
|
42 int i = mid - width;
|
Chris@2
|
43 if (i < 0)
|
Chris@2
|
44 i = 0;
|
Chris@2
|
45 int stop = mid + width + 1;
|
Chris@2
|
46 if (stop > data.length)
|
Chris@2
|
47 stop = data.length;
|
Chris@2
|
48 maxp = i;
|
Chris@2
|
49 for (i++; i < stop; i++)
|
Chris@2
|
50 if (data[i] > data[maxp])
|
Chris@2
|
51 maxp = i;
|
Chris@2
|
52 if (maxp == mid) {
|
Chris@2
|
53 int j;
|
Chris@2
|
54 for (j = peakCount; j > 0; j--) {
|
Chris@2
|
55 if (data[maxp] <= data[peaks[j-1]])
|
Chris@2
|
56 break;
|
Chris@2
|
57 else if (j < peaks.length)
|
Chris@2
|
58 peaks[j] = peaks[j-1];
|
Chris@2
|
59 }
|
Chris@2
|
60 if (j != peaks.length)
|
Chris@2
|
61 peaks[j] = maxp;
|
Chris@2
|
62 if (peakCount != peaks.length)
|
Chris@2
|
63 peakCount++;
|
Chris@2
|
64 }
|
Chris@2
|
65 mid++;
|
Chris@2
|
66 }
|
Chris@2
|
67 return peakCount;
|
Chris@2
|
68 } // findPeaks()
|
Chris@2
|
69
|
Chris@2
|
70 /** General peak picking method for finding local maxima in an array
|
Chris@2
|
71 * @param data input data
|
Chris@2
|
72 * @param width minimum distance between peaks
|
Chris@2
|
73 * @param threshold minimum value of peaks
|
Chris@2
|
74 * @return list of peak indexes
|
Chris@2
|
75 */
|
Chris@2
|
76 public static LinkedList<Integer> findPeaks(double[] data, int width,
|
Chris@2
|
77 double threshold) {
|
Chris@2
|
78 return findPeaks(data, width, threshold, 0, false);
|
Chris@2
|
79 } // findPeaks()
|
Chris@2
|
80
|
Chris@2
|
81 /** General peak picking method for finding local maxima in an array
|
Chris@2
|
82 * @param data input data
|
Chris@2
|
83 * @param width minimum distance between peaks
|
Chris@2
|
84 * @param threshold minimum value of peaks
|
Chris@2
|
85 * @param decayRate how quickly previous peaks are forgotten
|
Chris@2
|
86 * @param isRelative minimum value of peaks is relative to local average
|
Chris@2
|
87 * @return list of peak indexes
|
Chris@2
|
88 */
|
Chris@2
|
89 public static LinkedList<Integer> findPeaks(double[] data, int width,
|
Chris@2
|
90 double threshold, double decayRate, boolean isRelative) {
|
Chris@2
|
91 LinkedList<Integer> peaks = new LinkedList<Integer>();
|
Chris@2
|
92 int maxp = 0;
|
Chris@2
|
93 int mid = 0;
|
Chris@2
|
94 int end = data.length;
|
Chris@2
|
95 double av = data[0];
|
Chris@2
|
96 while (mid < end) {
|
Chris@2
|
97 av = decayRate * av + (1 - decayRate) * data[mid];
|
Chris@2
|
98 if (av < data[mid])
|
Chris@2
|
99 av = data[mid];
|
Chris@2
|
100 int i = mid - width;
|
Chris@2
|
101 if (i < 0)
|
Chris@2
|
102 i = 0;
|
Chris@2
|
103 int stop = mid + width + 1;
|
Chris@2
|
104 if (stop > data.length)
|
Chris@2
|
105 stop = data.length;
|
Chris@2
|
106 maxp = i;
|
Chris@2
|
107 for (i++; i < stop; i++)
|
Chris@2
|
108 if (data[i] > data[maxp])
|
Chris@2
|
109 maxp = i;
|
Chris@2
|
110 if (maxp == mid) {
|
Chris@2
|
111 if (overThreshold(data, maxp, width, threshold, isRelative,av)){
|
Chris@2
|
112 if (debug)
|
Chris@2
|
113 System.out.println(" peak");
|
Chris@2
|
114 peaks.add(new Integer(maxp));
|
Chris@2
|
115 } else if (debug)
|
Chris@2
|
116 System.out.println();
|
Chris@2
|
117 }
|
Chris@2
|
118 mid++;
|
Chris@2
|
119 }
|
Chris@2
|
120 return peaks;
|
Chris@2
|
121 } // findPeaks()
|
Chris@2
|
122
|
Chris@2
|
123 public static double expDecayWithHold(double av, double decayRate,
|
Chris@2
|
124 double[] data, int start, int stop) {
|
Chris@2
|
125 while (start < stop) {
|
Chris@2
|
126 av = decayRate * av + (1 - decayRate) * data[start];
|
Chris@2
|
127 if (av < data[start])
|
Chris@2
|
128 av = data[start];
|
Chris@2
|
129 start++;
|
Chris@2
|
130 }
|
Chris@2
|
131 return av;
|
Chris@2
|
132 } // expDecayWithHold()
|
Chris@2
|
133
|
Chris@2
|
134 public static boolean overThreshold(double[] data, int index, int width,
|
Chris@2
|
135 double threshold, boolean isRelative,
|
Chris@2
|
136 double av) {
|
Chris@2
|
137 if (debug)
|
Chris@2
|
138 System.out.printf("%4d : %6.3f Av1: %6.3f ",
|
Chris@2
|
139 index, data[index], av);
|
Chris@2
|
140 if (data[index] < av)
|
Chris@2
|
141 return false;
|
Chris@2
|
142 if (isRelative) {
|
Chris@2
|
143 int iStart = index - pre * width;
|
Chris@2
|
144 if (iStart < 0)
|
Chris@2
|
145 iStart = 0;
|
Chris@2
|
146 int iStop = index + post * width;
|
Chris@2
|
147 if (iStop > data.length)
|
Chris@2
|
148 iStop = data.length;
|
Chris@2
|
149 double sum = 0;
|
Chris@2
|
150 int count = iStop - iStart;
|
Chris@2
|
151 while (iStart < iStop)
|
Chris@2
|
152 sum += data[iStart++];
|
Chris@2
|
153 if (debug)
|
Chris@2
|
154 System.out.printf(" %6.3f %6.3f ", sum / count,
|
Chris@2
|
155 data[index] - sum / count - threshold);
|
Chris@2
|
156 return (data[index] > sum / count + threshold);
|
Chris@2
|
157 } else
|
Chris@2
|
158 return (data[index] > threshold);
|
Chris@2
|
159 } // overThreshold()
|
Chris@2
|
160
|
Chris@2
|
161 public static void normalise(double[] data) {
|
Chris@2
|
162 double sx = 0;
|
Chris@2
|
163 double sxx = 0;
|
Chris@2
|
164 for (int i = 0; i < data.length; i++) {
|
Chris@2
|
165 sx += data[i];
|
Chris@2
|
166 sxx += data[i] * data[i];
|
Chris@2
|
167 }
|
Chris@2
|
168 double mean = sx / data.length;
|
Chris@2
|
169 double sd = Math.sqrt((sxx - sx * mean) / data.length);
|
Chris@2
|
170 if (sd == 0)
|
Chris@2
|
171 sd = 1; // all data[i] == mean -> 0; avoids div by 0
|
Chris@2
|
172 for (int i = 0; i < data.length; i++) {
|
Chris@2
|
173 data[i] = (data[i] - mean) / sd;
|
Chris@2
|
174 }
|
Chris@2
|
175 } // normalise()
|
Chris@2
|
176
|
Chris@2
|
177 /** Uses an n-point linear regression to estimate the slope of data.
|
Chris@2
|
178 * @param data input data
|
Chris@2
|
179 * @param hop spacing of data points
|
Chris@2
|
180 * @param n length of linear regression
|
Chris@2
|
181 * @param slope output data
|
Chris@2
|
182 */
|
Chris@2
|
183 public static void getSlope(double[] data, double hop, int n,
|
Chris@2
|
184 double[] slope) {
|
Chris@2
|
185 int i = 0, j = 0;
|
Chris@2
|
186 double t;
|
Chris@2
|
187 double sx = 0, sxx = 0, sy = 0, sxy = 0;
|
Chris@2
|
188 for ( ; i < n; i++) {
|
Chris@2
|
189 t = i * hop;
|
Chris@2
|
190 sx += t;
|
Chris@2
|
191 sxx += t * t;
|
Chris@2
|
192 sy += data[i];
|
Chris@2
|
193 sxy += t * data[i];
|
Chris@2
|
194 }
|
Chris@2
|
195 double delta = n * sxx - sx * sx;
|
Chris@2
|
196 for ( ; j < n / 2; j++)
|
Chris@2
|
197 slope[j] = (n * sxy - sx * sy) / delta;
|
Chris@2
|
198 for ( ; j < data.length - (n + 1) / 2; j++, i++) {
|
Chris@2
|
199 slope[j] = (n * sxy - sx * sy) / delta;
|
Chris@2
|
200 sy += data[i] - data[i - n];
|
Chris@2
|
201 sxy += hop * (n * data[i] - sy);
|
Chris@2
|
202 }
|
Chris@2
|
203 for ( ; j < data.length; j++)
|
Chris@2
|
204 slope[j] = (n * sxy - sx * sy) / delta;
|
Chris@2
|
205 } // getSlope()
|
Chris@2
|
206
|
Chris@2
|
207 public static double min(double[] arr) { return arr[imin(arr)]; }
|
Chris@2
|
208
|
Chris@2
|
209 public static double max(double[] arr) { return arr[imax(arr)]; }
|
Chris@2
|
210
|
Chris@2
|
211 public static int imin(double[] arr) {
|
Chris@2
|
212 int i = 0;
|
Chris@2
|
213 for (int j = 1; j < arr.length; j++)
|
Chris@2
|
214 if (arr[j] < arr[i])
|
Chris@2
|
215 i = j;
|
Chris@2
|
216 return i;
|
Chris@2
|
217 } // imin()
|
Chris@2
|
218
|
Chris@2
|
219 public static int imax(double[] arr) {
|
Chris@2
|
220 int i = 0;
|
Chris@2
|
221 for (int j = 1; j < arr.length; j++)
|
Chris@2
|
222 if (arr[j] > arr[i])
|
Chris@2
|
223 i = j;
|
Chris@2
|
224 return i;
|
Chris@2
|
225 } // imax()
|
Chris@2
|
226
|
Chris@2
|
227 } // class Peaks
|