c@389
|
1 /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */
|
c@389
|
2
|
c@390
|
3 /*
|
c@390
|
4 QM DSP Library
|
c@390
|
5
|
c@390
|
6 Centre for Digital Music, Queen Mary, University of London.
|
c@390
|
7 This file Copyright 2010 Chris Cannam.
|
c@390
|
8
|
c@390
|
9 This program is free software; you can redistribute it and/or
|
c@390
|
10 modify it under the terms of the GNU General Public License as
|
c@390
|
11 published by the Free Software Foundation; either version 2 of the
|
c@390
|
12 License, or (at your option) any later version. See the file
|
c@390
|
13 COPYING included with this distribution for more information.
|
c@390
|
14 */
|
c@390
|
15
|
c@389
|
16 #ifndef MEDIAN_FILTER_H
|
c@389
|
17 #define MEDIAN_FILTER_H
|
c@389
|
18
|
c@389
|
19 #include <algorithm>
|
c@389
|
20 #include <cassert>
|
c@389
|
21 #include <cmath>
|
c@389
|
22 #include <iostream>
|
c@389
|
23
|
c@389
|
24 template <typename T>
|
c@389
|
25 class MedianFilter
|
c@389
|
26 {
|
c@389
|
27 public:
|
c@389
|
28 MedianFilter(int size, float percentile = 50.f) :
|
c@389
|
29 m_size(size),
|
c@389
|
30 m_frame(new T[size]),
|
c@389
|
31 m_sorted(new T[size]),
|
c@389
|
32 m_sortend(m_sorted + size - 1) {
|
c@389
|
33 setPercentile(percentile);
|
c@389
|
34 reset();
|
c@389
|
35 }
|
c@389
|
36
|
c@389
|
37 ~MedianFilter() {
|
c@389
|
38 delete[] m_frame;
|
c@389
|
39 delete[] m_sorted;
|
c@389
|
40 }
|
c@389
|
41
|
c@389
|
42 void setPercentile(float p) {
|
c@389
|
43 m_index = int((m_size * p) / 100.f);
|
c@389
|
44 if (m_index >= m_size) m_index = m_size-1;
|
c@389
|
45 if (m_index < 0) m_index = 0;
|
c@389
|
46 }
|
c@389
|
47
|
c@389
|
48 void push(T value) {
|
c@389
|
49 if (value != value) return; // nan
|
c@389
|
50 drop(m_frame[0]);
|
c@389
|
51 const int sz1 = m_size-1;
|
c@389
|
52 for (int i = 0; i < sz1; ++i) m_frame[i] = m_frame[i+1];
|
c@389
|
53 m_frame[m_size-1] = value;
|
c@389
|
54 put(value);
|
c@389
|
55 }
|
c@389
|
56
|
c@389
|
57 T get() const {
|
c@389
|
58 return m_sorted[m_index];
|
c@389
|
59 }
|
c@389
|
60
|
c@389
|
61 T getAt(float percentile) {
|
c@389
|
62 int ix = int((m_size * percentile) / 100.f);
|
c@389
|
63 if (ix >= m_size) ix = m_size-1;
|
c@389
|
64 if (ix < 0) ix = 0;
|
c@389
|
65 return m_sorted[ix];
|
c@389
|
66 }
|
c@389
|
67
|
c@389
|
68 void reset() {
|
c@389
|
69 for (int i = 0; i < m_size; ++i) m_frame[i] = 0;
|
c@389
|
70 for (int i = 0; i < m_size; ++i) m_sorted[i] = 0;
|
c@389
|
71 }
|
c@389
|
72
|
c@389
|
73 private:
|
c@389
|
74 const int m_size;
|
c@389
|
75 T *const m_frame;
|
c@389
|
76 T *const m_sorted;
|
c@389
|
77 T *const m_sortend;
|
c@389
|
78 int m_index;
|
c@389
|
79
|
c@389
|
80 void put(T value) {
|
c@389
|
81 // precondition: m_sorted contains m_size-1 values, packed at start
|
c@389
|
82 // postcondition: m_sorted contains m_size values, one of which is value
|
c@389
|
83 T *point = std::lower_bound(m_sorted, m_sortend, value);
|
c@389
|
84 const int n = m_sortend - point;
|
c@389
|
85 for (int i = n; i > 0; --i) point[i] = point[i-1];
|
c@389
|
86 *point = value;
|
c@389
|
87 }
|
c@389
|
88
|
c@389
|
89 void drop(T value) {
|
c@389
|
90 // precondition: m_sorted contains m_size values, one of which is value
|
c@389
|
91 // postcondition: m_sorted contains m_size-1 values, packed at start
|
c@389
|
92 T *point = std::lower_bound(m_sorted, m_sortend + 1, value);
|
c@389
|
93 if (*point != value) {
|
c@389
|
94 std::cerr << "WARNING: MedianFilter::drop: *point is " << *point
|
c@389
|
95 << ", expected " << value << std::endl;
|
c@389
|
96 }
|
c@389
|
97 const int n = m_sortend - point;
|
c@389
|
98 for (int i = 0; i < n; ++i) point[i] = point[i+1];
|
c@389
|
99 *m_sortend = T(0);
|
c@389
|
100 }
|
c@390
|
101
|
c@390
|
102 MedianFilter(const MedianFilter &); // not provided
|
c@390
|
103 MedianFilter &operator=(const MedianFilter &); // not provided
|
c@389
|
104 };
|
c@389
|
105
|
c@389
|
106 #endif
|
c@389
|
107
|