Recent

Author Topic: FFT scaling (old chestnut)  (Read 916 times)

MarkMLl

  • Hero Member
  • *****
  • Posts: 2088
FFT scaling (old chestnut)
« on: January 27, 2021, 11:59:31 am »
I wonder whether I could bore the community with a question which is more about using FPC/Lazarus effectively than the language or tools.

I have a cheap oscilloscope as a pocket tool: as described at https://hackaday.com/2017/08/08/touchscreen-oscilloscope/ it's a handy bit of kit for its price, although in practice it will need new leads with the appropriate (MMCX?) connectors.

I've written a program that captures data over USB, where- as is fairly common- the device looks like a serial port (CP2102). This is entirely adequate to save screenshots to go into project documentation etc., at least in the time domain. Running on Linux x86_64 with FPC 3.2.0 and Lazarus 2.1.0, but also tested on other versions/variants.

My problems start when I try to convert the collected time domain data to frequency domain using FFT. I was initially using the code discussed at e.g. https://lists.freepascal.org/pipermail/fpc-pascal/2017-April/050717.html but have also got the same questionable results using FFTW, so it's obviously my interpretation of the results that's at fault.

The scope's calibration signal is 1kHz 3.3V. If I capture that with the timebase set to 1 mSec per division I get 512 values, and can derive appropriate voltages etc. for on-screen annotation. There's 20 divisions, so it's 20 mSec across the entire 512 values or 39.06 uSec per value hence a sampling rate of 25.6 kHz.

/If/ my understanding is correct, that implies that each bin in the FFT output is 25.6 kHz/512, or 50Hz. However because this is mirrored with the zero in the middle and +-Nyquist at each end what I've actually got is a usable 256 bins from 0 to 12.8 kHz (still 50Hz per bin). Displaying that with an arbitrary 32 divisions should give me 400 Hz per division.

Please could somebody who knows what they're doing confirm my logic there before I wade into my code again. I attach a couple of screenshots, showing the sort of problem I'm having.

The bigger picture is that I'm trying to get the graphical side of this wrapped up as a way of revisiting TAChart. Then I'm going to add graphical output to a program that reads data from a multimeter (Mastech MS2115B), and the code that does that will be adapted in short order to read a pulse oximeter which is obviously of topical interest.

MarkMLl







Turbo Pascal v1 on CCP/M-86, multitasking with LAN and graphics in 128Kb.
Pet hate: people who boast about the size and sophistication of their computer.

Mr.Madguy

  • Hero Member
  • *****
  • Posts: 627
Re: FFT scaling (old chestnut)
« Reply #1 on: January 27, 2021, 12:52:54 pm »
I'm not sure about math, because I'm not strong in it, but I use following code:
Code: Pascal  [Select][+][-]
  1. unit FFT;
  2.  
  3. {$mode objfpc}{$H+}
  4.  
  5. interface
  6.  
  7. uses Complex;
  8.  
  9. procedure FFTMain(ADepth:Integer;ASource, ADest:PComplexArray);
  10.  
  11. implementation
  12.  
  13. procedure DoFFT(ADepth, AMaxDepth:Integer;ASource, ADest:PComplexArray;ATable:PComplexArray);
  14.   var Shift, SourceSpacing, DestSpacing, I:Integer;Temp:TComplex;
  15. begin
  16.   if ADepth = 0 then begin
  17.     ADest^[0] := ASource^[0];
  18.     Exit;
  19.   end;
  20.   Shift := AMaxDepth - ADepth;
  21.   SourceSpacing := 1 shl Shift;
  22.   DestSpacing := 1 shl (ADepth - 1);
  23.   DoFFT(ADepth - 1, AMaxDepth, ASource, ADest, ATable);
  24.   DoFFT(ADepth - 1, AMaxDepth, @ASource^[SourceSpacing], @ADest^[DestSpacing], ATable);
  25.   for I := 0 to DestSpacing - 1 do begin
  26.     Temp := ADest^[DestSpacing + I] * ATable^[I shl Shift];
  27.     ADest^[DestSpacing + I] := ADest^[I] - Temp;
  28.     ADest^[I] := ADest^[I] + Temp;
  29.   end;
  30. end;
  31.  
  32. procedure FFTMain(ADepth:Integer;ASource, ADest:PComplexArray);
  33.   var Table:PComplexArray;Count, I:Integer;
  34. begin
  35.   Count := 1 shl ADepth;
  36.   Table := GetMem(Count * SizeOf(TComplex));
  37.   for I := 0 to Count - 1 do begin
  38.     Table^[I].Re := Cos(-(2 * Pi) * I / Count);
  39.     Table^[I].Im := Sin(-(2 * Pi) * I / Count);
  40.   end;
  41.   DoFFT(ADepth, ADepth, ASource, ADest, Table);
  42.   for I := 0 to Count - 1 do begin
  43.     ADest^[I] := ADest^[I] / Count;
  44.   end;
  45.   FreeMem(Table);
  46. end;
  47.  
  48. end.
  49.  
And only problem I have - is mirroring. Picture 1 shows "left" mirroring and picture 2 shows "right" one. I don't know, how to beat it and I don't care, because it's just example of uncertainty principle in quantum mechanic.
« Last Edit: January 27, 2021, 12:57:37 pm by Mr.Madguy »
DynamicData 3.0 is released!
4.0 version is in development! Another complete overhaul! Things are made more simple and, of course, faster!
Aiming at Lazarus support, when closures will be implemented.

avra

  • Hero Member
  • *****
  • Posts: 2112
    • Additional info
Re: FFT scaling (old chestnut)
« Reply #2 on: January 27, 2021, 01:15:02 pm »
...hence a sampling rate of 25.6 kHz.
This does not sound right. Can you test with 0.1 ms/div and report difference in findings?

How do you read data? Do you get CSV file as explained at http://www.jyetech.com/Products/LcdScope/CaptureUploadingAndWaveformFileFormat.pdf or you use USB Scope mode and DataBlock read function as in http://www.jyetech.com/Products/LcdScope/TheDataInterface_068.pdf ? Do you use Auto state or Manual state?

You can compare FFT results of your sample data online https://duckduckgo.com/?q=fft+online&ia=web.

Btw, there is an unofficial sigrok driver for your device:
https://sourceforge.net/p/sigrok/mailman/message/35238342/
https://github.com/mwm/libsigrok/tree/dso112a/src/hardware/jyetech-dso112a

EDIT:

From Data Interface manual:
Quote
Note that when the timebase is set between 20ms/div – 0.5us/div (inclusive) the whole buffer
of data is sent (as DataBlock frame) after each capture. When the timebase is set between
10min/div and 50ms/div (inclusive) single sample is sent (as DataSample frame) after each
sampling.
Maybe buffer length is 512 and your perception is wrong? If you get the same 512 bytes as block size for both 1 ms/div and 0.1 ms/div then it might be the case. Or device always sends fixed 512 bytes for 20 divisions?
« Last Edit: January 27, 2021, 01:52:51 pm by avra »
ct2laz - Conversion between Lazarus and CodeTyphon
bithelpers - Bit manipulation for standard types
pasettimino - Siemens S7 PLC lib

MarkMLl

  • Hero Member
  • *****
  • Posts: 2088
Re: FFT scaling (old chestnut)
« Reply #3 on: January 27, 2021, 01:28:49 pm »
...hence a sampling rate of 25.6 kHz.
This does not sound right. Can you test with 0.1 ms/div and report difference in findings?

That was a worked-out figure rather than an observed one, so I'm very interested to know /why/ it doesn't sound right since if I've got it wrong there it might explain what's wrong.

Quote
How do you read data? Do you get CSV file as explained at http://www.jyetech.com/Products/LcdScope/CaptureUploadingAndWaveformFileFormat.pdf or you use USB Scope mode and DataBlock read function as in http://www.jyetech.com/Products/LcdScope/TheDataInterface_068.pdf ?

CSV via Xmodem, since other devices I hope to support at some point only have a "dump" rather than a "live-control" facility and much of this is already in a separate dynamically-linked protocol unit.

Quote
You can compare FFT results of your sample data online https://duckduckgo.com/?q=fft+online&ia=web.

Thanks. I'll check.

Quote
Btw, there is an unofficial sigrok driver for your device:
https://sourceforge.net/p/sigrok/mailman/message/35238342/
https://github.com/mwm/libsigrok/tree/dso112a/src/hardware/jyetech-dso112a

Yes, but some of the other devices I hope to support aren't so I thought it worth using this as (in part) a testbed for dynamic linkage stuff I've had running in the past in a different context. Problems raised by that aren't exactly on the critical path at the moment, so I'll raise them some other time :-)

MarkMLl
Turbo Pascal v1 on CCP/M-86, multitasking with LAN and graphics in 128Kb.
Pet hate: people who boast about the size and sophistication of their computer.

avra

  • Hero Member
  • *****
  • Posts: 2112
    • Additional info
Re: FFT scaling (old chestnut)
« Reply #4 on: January 27, 2021, 02:05:16 pm »
That was a worked-out figure rather than an observed one, so I'm very interested to know /why/ it doesn't sound right since if I've got it wrong there it might explain what's wrong.

I have seen somewhere that they used AVR MEGA64 in some of their DSOs, so I will assume that your model uses AVR too. 2.5MHz sample rate is not achievable with AVR internal ADC, so I will also assume that they use some fast external ADC. That ADC must be tied to some clock, and 25.6 kHz could be achieved from something like 16MHz clock by using some prescaler. That could be partially determined by opening device and examination of board components and ADC datasheet, but I do not have a device to examine so I am just guessing. I have also updated my first post with a question if datablock has fixed size of 512 bytes but you have replied before my edit was finished so you probably missed it. You could also check if you can try to receive CSV with TeraTerm as described in CaptureUploadingAndWaveformFileFormat.pdf because that CSV shows RecLen and SampleRate for the transferred data which could be used for clarification of your doubt.

Btw. SampleRate in that pdf is shown to be 50000, so if that is right then your 25600 is questionable. Also RecLen is 256, so fixed buffer of 512 bytes is also questionable. Again, only if that DSO068 screenshot is real and applicable to your DSO model (I have seen some forum message which says that DSO068 protocol  is compatible with your DSO).
« Last Edit: January 27, 2021, 02:08:48 pm by avra »
ct2laz - Conversion between Lazarus and CodeTyphon
bithelpers - Bit manipulation for standard types
pasettimino - Siemens S7 PLC lib

wp

  • Hero Member
  • *****
  • Posts: 8107
Re: FFT scaling (old chestnut)
« Reply #5 on: January 27, 2021, 03:20:54 pm »
The bigger picture is that I'm trying to get the graphical side of this wrapped up as a way of revisiting TAChart. Then I'm going to add graphical output to a program that reads data from a multimeter (Mastech MS2115B), and the code that does that will be adapted in short order to read a pulse oximeter which is obviously of topical interest.
Your projects seem to be similar to two of mine: oscilloscope reading data from the audio input connector (using the BASS libarary; https://sourceforge.net/p/wp-laz/code/HEAD/tree/Oscilloscope/), and logging analog values from a digital multimeter based on the devices FS9721 and FS9922 (https://sourceforge.net/p/wp-laz/code/HEAD/tree/datalogger/). I am using TAChart a lot, and maybe you can find some ideas there.
Mainly Lazarus trunk / fpc 3.2.0 / all 32-bit on Win-10, but many more...

MarkMLl

  • Hero Member
  • *****
  • Posts: 2088
Re: FFT scaling (old chestnut)
« Reply #6 on: January 27, 2021, 03:42:04 pm »
That was a worked-out figure rather than an observed one, so I'm very interested to know /why/ it doesn't sound right since if I've got it wrong there it might explain what's wrong.

I have seen somewhere that they used AVR MEGA64 in some of their DSOs, so I will assume that your model uses AVR too. 2.5MHz sample rate is not achievable with AVR internal ADC, so I will also assume that they use some fast external ADC. That ADC must be tied to some clock, and 25.6 kHz could be achieved from something like 16MHz clock by using some prescaler. That could be partially determined by opening device and examination of board components and ADC datasheet, but I do not have a device to examine so I am just guessing.

Photo at http://www.hezkyden.cz/i/osciloskop-dso112a-konstrukce.jpg plus I've got mine open. 2x Atmels (I've previously remarked on this since the loader is written using Delphi), and various other stuff under a non-removable screening PCB. There's a 20MHz xtal in there but I suspect that that's on an external circuit supplying both chips... HP bench 'scope is happy enough with that being a 20MHz signal.

Under the circumstances I'm not bothered with it being out by a few percent, just trying to get it better than a factor of 2.
[/quote]

Quote
I have also updated my first post with a question if datablock has fixed size of 512 bytes but you have replied before my edit was finished so you probably missed it. You could also check if you can try to receive CSV with TeraTerm as described in CaptureUploadingAndWaveformFileFormat.pdf because that CSV shows RecLen and SampleRate for the transferred data which could be used for clarification of your doubt.

Who is this Teraterm of which you speak? :-) OK, I know perfectly well what it is but I'm running Linux as stated in the first message. Download debug messages confirm 512 readings (not bytes) being transferred, but I've found a slightly dodgy assumption relating to the time-per-reading that I need to check.

Quote
Btw. SampleRate in that pdf is shown to be 50000, so if that is right then your 25600 is questionable. Also RecLen is 256, so fixed buffer of 512 bytes is also questionable. Again, only if that DSO068 screenshot is real and applicable to your DSO model (I have seen some forum message which says that DSO068 protocol  is compatible with your DSO).

But it's XModem and transferring multiple records... again I can see the entire session in debug messages:


...
> 06
< 01 1A E5
< 30 38 38 0A 30 30 30 38 38 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 30 38 38 0A 30 30 30 38 37 0A 30 30 30 38 38 0A 30 30 30 38 37 0A 30 30 30 38 38 0A 30 30 30 38 37 0A 30 30 30 38 38 0A 30 30 30 38
< D7
> 06
< 01 1B E4
< 37 0A 30 30 30 38 38 0A 30 30 30 38 37 0A 30 30 30 38 38 0A 30 30 30 38 37 0A 30 30 30 38 38 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 31 37 30 0A 30 30 30 38 37 0A 30 30 30 38 38 0A 30 30 30 38 37 0A 30 30 30 38 38 0A
< AF
> 06
< 04 00 00
> 06


Note also that I think that that edit you put into your original message applies to live data, not for a dump (which needs action both on the PC and on the 'scope, and is a one-shot operation).

So leaving aside for the moment that my understanding of the sample rate might be a few percent adrift, do you see any howlers in my original logic? If it looks more or less OK I'll fix my frequency domain bin size to conform to it, and will then review my assumptions about time domain scaling... I suspect that I'm doing something like neglecting the size of the display window when working out the samples-per-division figure.

MarkMLl
Turbo Pascal v1 on CCP/M-86, multitasking with LAN and graphics in 128Kb.
Pet hate: people who boast about the size and sophistication of their computer.

avra

  • Hero Member
  • *****
  • Posts: 2112
    • Additional info
Re: FFT scaling (old chestnut)
« Reply #7 on: January 27, 2021, 03:55:18 pm »
16MHz can be divided with 25.6 kHz, but 20MHz can not, so if you can not find any other crystals except already mentioned 20MHz, I highly doubt that sampling rate is 25600.

As mentioned, CaptureUploadingAndWaveformFileFormat.pdf on page 3 figure 7 (which I have attached here) shows CSV received with TeraTerrm and procedure to get is described. I know it's windows only, but I do not see any other way to determine exact SampleRate and RecLen filelds for received data. Once you get them, you can compare to your assumptions and make corrections.
ct2laz - Conversion between Lazarus and CodeTyphon
bithelpers - Bit manipulation for standard types
pasettimino - Siemens S7 PLC lib

MarkMLl

  • Hero Member
  • *****
  • Posts: 2088
Re: FFT scaling (old chestnut)
« Reply #8 on: January 27, 2021, 04:06:52 pm »
The bigger picture is that I'm trying to get the graphical side of this wrapped up as a way of revisiting TAChart. Then I'm going to add graphical output to a program that reads data from a multimeter (Mastech MS2115B), and the code that does that will be adapted in short order to read a pulse oximeter which is obviously of topical interest.
Your projects seem to be similar to two of mine: oscilloscope reading data from the audio input connector (using the BASS libarary; https://sourceforge.net/p/wp-laz/code/HEAD/tree/Oscilloscope/), and logging analog values from a digital multimeter based on the devices FS9721 and FS9922 (https://sourceforge.net/p/wp-laz/code/HEAD/tree/datalogger/). I am using TAChart a lot, and maybe you can find some ideas there.

Thanks, I noticed your comment about TAChart a few days ago and intended to follow up your link. In my case I've not so much got end-user applications in mind but am trying to wring every bit of functionality out of what-I've got test equipment... hence my occasional comment about GPIB etc. :-)

I can't remember what chip's in the 2115B although I'm obviously aware that there's a limited number used in practice. I chose that meter because it had the useful combination of a current clamp and data download capability (albeit with a fairly slow sample rate); data format is a bit more advanced than what's described in your PDFs with lots of info on scaling etc.

I uploaded a text-only program to read it to Github over the weekend (and then had to correct it due to the FileExists() issue which was messing up the code that worked out what port it was on) https://github.com/MarkMLl/Mastech_ms2115b In practical terms very few changes will be needed to handle a cheap pulse oximeter.

MarkMLl
Turbo Pascal v1 on CCP/M-86, multitasking with LAN and graphics in 128Kb.
Pet hate: people who boast about the size and sophistication of their computer.

MarkMLl

  • Hero Member
  • *****
  • Posts: 2088
Re: FFT scaling (old chestnut)
« Reply #9 on: January 27, 2021, 04:09:15 pm »
16MHz can be divided with 25.6 kHz, but 20MHz can not, so if you can not find any other crystals except already mentioned 20MHz, I highly doubt that sampling rate is 25600.

As mentioned, CaptureUploadingAndWaveformFileFormat.pdf on page 3 figure 7 (which I have attached here) shows CSV received with TeraTerrm and procedure to get is described. I know it's windows only, but I do not see any other way to determine exact SampleRate and RecLen filelds for received data. Once you get them, you can compare to your assumptions and make corrections.

Remember that that's not the same kind of meter so the fact that the number of samples differs is hardly surprising (it's also got a 1024-sample mode but I've not tried that at all with FFT etc. yet). I've got no problem at all capturing the data, the problem's in my interpretation when FFTing it.

MarkMLl
« Last Edit: January 27, 2021, 04:24:23 pm by MarkMLl »
Turbo Pascal v1 on CCP/M-86, multitasking with LAN and graphics in 128Kb.
Pet hate: people who boast about the size and sophistication of their computer.

Fred vS

  • Hero Member
  • *****
  • Posts: 2085
    • StrumPract is the musicians best friend
Re: FFT scaling (old chestnut)
« Reply #10 on: January 27, 2021, 05:52:34 pm »
Your projects seem to be similar to two of mine: oscilloscope reading data from the audio input connector (using the BASS libarary; https://sourceforge.net/p/wp-laz/code/HEAD/tree/Oscilloscope/)

Hello Wp.

Nice project!

I try to compile it but there is a missing component: uEControls.
Where can I find the component?

Your project seems to be Windows only.
What would it needed to make it work under Unix OS (like Raspbian for Rpi) ?

Your project uses Bass closed source library.
What methods of Bass are you using, do you think that open source libraries like PortAudio + LibSndFile could do the job?

Thanks.

Fre;D
I use Lazarus 2.0.6 32/64 and FPC 3.2.0 32/64 on Debian 10.2 64 bit, Windows 10, Windows 7 32/64, Windows XP 32,  FreeBSD 64 and Mac OS X Snow Leopard 32.
Widgetset: fpGUI, MSEgui, Win32, GTK2, Qt, Carbon.

https://github.com/fredvs
https://gitlab.com/fredvs

wp

  • Hero Member
  • *****
  • Posts: 8107
Re: FFT scaling (old chestnut)
« Reply #11 on: January 27, 2021, 06:07:55 pm »
"Oscilloscope" was some kind of fun-project just to show the feasiblity. Don't expect too much of it.

UEControls is available via OPM.

Windows-only? Really? Where is it not cross-platform? It should be, although I never tested it.

Other audio libs will certainly work as well. I selected BASS because I had used it in an old Delphi project and was trying to avoid learning how to install and apply other libs. The "Oscilloscope" application uses some kind of abstract TDataCollector class to access the sound card. For BASS a special TBASSDataCollector is implemented. Likewise any another audio lib could be used provided  that a TDatacollector descendant can be written for it. Again: not severely tested, probably needing some tweaks here and there.
Mainly Lazarus trunk / fpc 3.2.0 / all 32-bit on Win-10, but many more...

Fred vS

  • Hero Member
  • *****
  • Posts: 2085
    • StrumPract is the musicians best friend
Re: FFT scaling (old chestnut)
« Reply #12 on: January 27, 2021, 06:34:34 pm »
UEControls is available via OPM.

Ha, ok, thanks.

Windows-only? Really? Where is it not cross-platform? It should be, although I never tested it.

So, no "Windows only" things, ok, perfect.

The "Oscilloscope" application uses some kind of abstract TDataCollector class

OK, I will re-jump into Bass code to remember what it is about.

Many thanks.

PS: No plan to migrate to GitHub (easier for pull-request) ?
I use Lazarus 2.0.6 32/64 and FPC 3.2.0 32/64 on Debian 10.2 64 bit, Windows 10, Windows 7 32/64, Windows XP 32,  FreeBSD 64 and Mac OS X Snow Leopard 32.
Widgetset: fpGUI, MSEgui, Win32, GTK2, Qt, Carbon.

https://github.com/fredvs
https://gitlab.com/fredvs

MarkMLl

  • Hero Member
  • *****
  • Posts: 2088
Re: FFT scaling (old chestnut)
« Reply #13 on: January 27, 2021, 07:11:58 pm »
I think I've found out what's going on. In the time domain, I was displaying a fixed number of horizontal divisions (20) and saying that the the timing (secs per division) was the same as that reported by the 'scope. However, I was also displaying the entire buffer: 512 or 1024 values rather than the 25x12=300 (I think) values in the 'scope's display window: if I put a label on the waveform the time and voltage were both plausible, but the horizontal spacing of the grid was slightly out when working with 512 values and grossly out with 1024.

Because of the way I was working with divisions, when I started converting from time divisions to frequency bins things went badly wrong.

All of that is obviously fixable, but it means I have to go back and get my time domain grid right before I can really tackle the frequency domain.

In order to keep things looking like a traditional scope, I'm seriously considering dropping the rightmost few readings, since otherwise I'll either have to move the timing divisions away from the conventional 1.000 mSec settings or won't have a grid line in the middle of the display as is traditional.

MarkMLl
Turbo Pascal v1 on CCP/M-86, multitasking with LAN and graphics in 128Kb.
Pet hate: people who boast about the size and sophistication of their computer.

wp

  • Hero Member
  • *****
  • Posts: 8107
Re: FFT scaling (old chestnut)
« Reply #14 on: January 27, 2021, 07:33:58 pm »
No plan to migrate to GitHub (easier for pull-request) ?
Here you go: https://github.com/wp-xyz/Oscilloscope.git
Mainly Lazarus trunk / fpc 3.2.0 / all 32-bit on Win-10, but many more...

 

TinyPortal © 2005-2018