1 /*
2  *  Copyright (C) 2015, Mike Walters <mike@flomp.net>
3  *
4  *  This file is part of inspectrum.
5  *
6  *  This program is free software: you can redistribute it and/or modify
7  *  it under the terms of the GNU General Public License as published by
8  *  the Free Software Foundation, either version 3 of the License, or
9  *  (at your option) any later version.
10  *
11  *  This program is distributed in the hope that it will be useful,
12  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
14  *  GNU General Public License for more details.
15  *
16  *  You should have received a copy of the GNU General Public License
17  *  along with this program.  If not, see <http://www.gnu.org/licenses/>.
18  */
19 
20 #include "spectrogramplot.h"
21 
22 #include <QDebug>
23 #include <QElapsedTimer>
24 #include <QPainter>
25 #include <QPaintEvent>
26 #include <QPixmapCache>
27 #include <QRect>
28 #include <liquid/liquid.h>
29 #include <functional>
30 #include <cstdlib>
31 #include "util.h"
32 
33 
SpectrogramPlot(std::shared_ptr<SampleSource<std::complex<float>>> src)34 SpectrogramPlot::SpectrogramPlot(std::shared_ptr<SampleSource<std::complex<float>>> src) : Plot(src), inputSource(src), fftSize(512), tuner(fftSize, this)
35 {
36     setFFTSize(fftSize);
37     zoomLevel = 1;
38     powerMax = 0.0f;
39     powerMin = -50.0f;
40     sampleRate = 0;
41     frequencyScaleEnabled = false;
42 
43     for (int i = 0; i < 256; i++) {
44         float p = (float)i / 256;
45         colormap[i] = QColor::fromHsvF(p * 0.83f, 1.0, 1.0 - p).rgba();
46     }
47 
48     tunerTransform = std::make_shared<TunerTransform>(src);
49     connect(&tuner, &Tuner::tunerMoved, this, &SpectrogramPlot::tunerMoved);
50 }
51 
invalidateEvent()52 void SpectrogramPlot::invalidateEvent()
53 {
54     pixmapCache.clear();
55     fftCache.clear();
56     emit repaint();
57 }
58 
paintFront(QPainter & painter,QRect & rect,range_t<size_t> sampleRange)59 void SpectrogramPlot::paintFront(QPainter &painter, QRect &rect, range_t<size_t> sampleRange)
60 {
61     if (tunerEnabled())
62         tuner.paintFront(painter, rect, sampleRange);
63 
64     if (frequencyScaleEnabled)
65         paintFrequencyScale(painter, rect);
66 }
67 
paintFrequencyScale(QPainter & painter,QRect & rect)68 void SpectrogramPlot::paintFrequencyScale(QPainter &painter, QRect &rect)
69 {
70     // At which pixel is F_+sampleRate/2
71     int y = rect.y();
72 
73     int plotHeight = rect.height();
74 
75     double bwPerPixel = (double)sampleRate / plotHeight;
76     int tickHeight = 50;
77 
78     int bwPerTick = 10 * pow(10, floor(log(bwPerPixel * tickHeight) / log(10)));
79 
80     if (bwPerTick < 1) {
81         return;
82     }
83 
84     painter.save();
85 
86     QPen pen(Qt::white, 1, Qt::SolidLine);
87     painter.setPen(pen);
88     QFontMetrics fm(painter.font());
89 
90 
91     int tick = 0;
92 
93     while (tick <= sampleRate / 2) {
94 
95         int tickpy = plotHeight / 2 - tick / bwPerPixel + y;
96         int tickny = plotHeight / 2 + tick / bwPerPixel + y;
97 
98         painter.drawLine(0, tickny, 30, tickny);
99         painter.drawLine(0, tickpy, 30, tickpy);
100 
101         if (tick != 0) {
102             char buf[128];
103 
104             if (bwPerTick % 1000000 == 0) {
105                 snprintf(buf, sizeof(buf), "-%d MHz", (int)tick / 1000000);
106             } else if(bwPerTick % 1000 == 0) {
107                 snprintf(buf, sizeof(buf), "-%d kHz", tick / 1000);
108             } else {
109                 snprintf(buf, sizeof(buf), "-%d Hz", tick);
110             }
111 
112             painter.drawText(5, tickny - 5, buf);
113 
114             buf[0] = ' ';
115             painter.drawText(5, tickpy + 15, buf);
116         }
117 
118         tick += bwPerTick;
119     }
120 
121     // Draw small ticks
122     bwPerTick /= 10;
123 
124     if (bwPerTick >= 1 ) {
125         tick = 0;
126         while (tick <= sampleRate / 2) {
127 
128             int tickpy = plotHeight / 2 - tick / bwPerPixel + y;
129             int tickny = plotHeight / 2 + tick / bwPerPixel + y;
130 
131             painter.drawLine(0, tickny, 3, tickny);
132             painter.drawLine(0, tickpy, 3, tickpy);
133 
134             tick += bwPerTick;
135         }
136     }
137     painter.restore();
138 }
139 
paintMid(QPainter & painter,QRect & rect,range_t<size_t> sampleRange)140 void SpectrogramPlot::paintMid(QPainter &painter, QRect &rect, range_t<size_t> sampleRange)
141 {
142     if (!inputSource || inputSource->count() == 0)
143         return;
144 
145     size_t sampleOffset = sampleRange.minimum % (getStride() * linesPerTile());
146     size_t tileID = sampleRange.minimum - sampleOffset;
147     int xoffset = sampleOffset / getStride();
148 
149     // Paint first (possibly partial) tile
150     painter.drawPixmap(QRect(rect.left(), rect.y(), linesPerTile() - xoffset, fftSize), *getPixmapTile(tileID), QRect(xoffset, 0, linesPerTile() - xoffset, fftSize));
151     tileID += getStride() * linesPerTile();
152 
153     // Paint remaining tiles
154     for (int x = linesPerTile() - xoffset; x < rect.right(); x += linesPerTile()) {
155         // TODO: don't draw past rect.right()
156         // TODO: handle partial final tile
157         painter.drawPixmap(QRect(x, rect.y(), linesPerTile(), fftSize), *getPixmapTile(tileID));
158         tileID += getStride() * linesPerTile();
159     }
160 }
161 
getPixmapTile(size_t tile)162 QPixmap* SpectrogramPlot::getPixmapTile(size_t tile)
163 {
164     QPixmap *obj = pixmapCache.object(TileCacheKey(fftSize, zoomLevel, tile));
165     if (obj != 0)
166         return obj;
167 
168     float *fftTile = getFFTTile(tile);
169     obj = new QPixmap(linesPerTile(), fftSize);
170     QImage image(linesPerTile(), fftSize, QImage::Format_RGB32);
171     float powerRange = -1.0f / std::abs(int(powerMin - powerMax));
172     for (int y = 0; y < fftSize; y++) {
173         auto scanLine = (QRgb*)image.scanLine(fftSize - y - 1);
174         for (int x = 0; x < linesPerTile(); x++) {
175             float *fftLine = &fftTile[x * fftSize];
176             float normPower = (fftLine[y] - powerMax) * powerRange;
177             normPower = clamp(normPower, 0.0f, 1.0f);
178 
179             scanLine[x] = colormap[(uint8_t)(normPower * (256 - 1))];
180         }
181     }
182     obj->convertFromImage(image);
183     pixmapCache.insert(TileCacheKey(fftSize, zoomLevel, tile), obj);
184     return obj;
185 }
186 
getFFTTile(size_t tile)187 float* SpectrogramPlot::getFFTTile(size_t tile)
188 {
189     std::array<float, tileSize>* obj = fftCache.object(TileCacheKey(fftSize, zoomLevel, tile));
190     if (obj != nullptr)
191         return obj->data();
192 
193     std::array<float, tileSize>* destStorage = new std::array<float, tileSize>;
194     float *ptr = destStorage->data();
195     size_t sample = tile;
196     while ((ptr - destStorage->data()) < tileSize) {
197         getLine(ptr, sample);
198         sample += getStride();
199         ptr += fftSize;
200     }
201     fftCache.insert(TileCacheKey(fftSize, zoomLevel, tile), destStorage);
202     return destStorage->data();
203 }
204 
getLine(float * dest,size_t sample)205 void SpectrogramPlot::getLine(float *dest, size_t sample)
206 {
207     if (inputSource && fft) {
208         auto buffer = inputSource->getSamples(sample, fftSize);
209         if (buffer == nullptr) {
210             auto neg_infinity = -1 * std::numeric_limits<float>::infinity();
211             for (int i = 0; i < fftSize; i++, dest++)
212                 *dest = neg_infinity;
213             return;
214         }
215 
216         for (int i = 0; i < fftSize; i++) {
217             buffer[i] *= window[i];
218         }
219 
220         fft->process(buffer.get(), buffer.get());
221         const float invFFTSize = 1.0f / fftSize;
222         const float logMultiplier = 10.0f / log2f(10.0f);
223         for (int i = 0; i < fftSize; i++) {
224             // Start from the middle of the FFTW array and wrap
225             // to rearrange the data
226             int k = i ^ (fftSize >> 1);
227             auto s = buffer[k] * invFFTSize;
228             float power = s.real() * s.real() + s.imag() * s.imag();
229             float logPower = log2f(power) * logMultiplier;
230             *dest = logPower;
231             dest++;
232         }
233     }
234 }
235 
getStride()236 int SpectrogramPlot::getStride()
237 {
238     return fftSize / zoomLevel;
239 }
240 
getTunerPhaseInc()241 float SpectrogramPlot::getTunerPhaseInc()
242 {
243     auto freq = 0.5f - tuner.centre() / (float)fftSize;
244     return freq * Tau;
245 }
246 
getTunerTaps()247 std::vector<float> SpectrogramPlot::getTunerTaps()
248 {
249     float cutoff = tuner.deviation() / (float)fftSize;
250     float gain = pow(10.0f, powerMax / -10.0f);
251     auto atten = 60.0f;
252     auto len = estimate_req_filter_len(0.05f, atten);
253     auto taps = std::vector<float>(len);
254     liquid_firdes_kaiser(len, cutoff, atten, 0.0f, taps.data());
255     std::transform(taps.begin(), taps.end(), taps.begin(),
256                    std::bind1st(std::multiplies<float>(), gain));
257     return taps;
258 }
259 
linesPerTile()260 int SpectrogramPlot::linesPerTile()
261 {
262     return tileSize / fftSize;
263 }
264 
mouseEvent(QEvent::Type type,QMouseEvent event)265 bool SpectrogramPlot::mouseEvent(QEvent::Type type, QMouseEvent event)
266 {
267     if (tunerEnabled())
268         return tuner.mouseEvent(type, event);
269 
270     return false;
271 }
272 
output()273 std::shared_ptr<AbstractSampleSource> SpectrogramPlot::output()
274 {
275     return tunerTransform;
276 }
277 
setFFTSize(int size)278 void SpectrogramPlot::setFFTSize(int size)
279 {
280     float sizeScale = float(size) / float(fftSize);
281     fftSize = size;
282     fft.reset(new FFT(fftSize));
283 
284     window.reset(new float[fftSize]);
285     for (int i = 0; i < fftSize; i++) {
286         window[i] = 0.5f * (1.0f - cos(Tau * i / (fftSize - 1)));
287     }
288 
289     setHeight(fftSize);
290     auto dev = tuner.deviation();
291     auto centre = tuner.centre();
292     tuner.setHeight(fftSize);
293     tuner.setDeviation( dev * sizeScale );
294     tuner.setCentre( centre * sizeScale );
295 }
296 
setPowerMax(int power)297 void SpectrogramPlot::setPowerMax(int power)
298 {
299     powerMax = power;
300     pixmapCache.clear();
301     tunerMoved();
302 }
303 
setPowerMin(int power)304 void SpectrogramPlot::setPowerMin(int power)
305 {
306     powerMin = power;
307     pixmapCache.clear();
308 }
309 
setZoomLevel(int zoom)310 void SpectrogramPlot::setZoomLevel(int zoom)
311 {
312     zoomLevel = zoom;
313 }
314 
setSampleRate(size_t rate)315 void SpectrogramPlot::setSampleRate(size_t rate)
316 {
317     sampleRate = rate;
318 }
319 
enableScales(bool enabled)320 void SpectrogramPlot::enableScales(bool enabled)
321 {
322    frequencyScaleEnabled = enabled;
323 }
324 
tunerEnabled()325 bool SpectrogramPlot::tunerEnabled()
326 {
327     return (tunerTransform->subscriberCount() > 0);
328 }
329 
tunerMoved()330 void SpectrogramPlot::tunerMoved()
331 {
332     tunerTransform->setFrequency(getTunerPhaseInc());
333     tunerTransform->setTaps(getTunerTaps());
334     tunerTransform->setRelativeBandwith(tuner.deviation() * 2.0 / getStride());
335 
336     // TODO: for invalidating traceplot cache, this shouldn't really go here
337     QPixmapCache::clear();
338 
339     emit repaint();
340 }
341 
qHash(const TileCacheKey & key,uint seed)342 uint qHash(const TileCacheKey &key, uint seed)
343 {
344     return key.fftSize ^ key.zoomLevel ^ key.sample ^ seed;
345 }
346