in reply to Re: Quick script to check data logger data
in thread Quick script to check data logger data

I pick PDL up for tasks occasionally, but not often enough to maintain fluency with it. It's one of those catch-22 situations where PDL would be a really handy tool to have in the toolbox, but <fluency> * <need to use> never quite peeks over the <is more productive with> threshold. Just a little more <need to use> would increase <fluency> enough to make the difference, but I find too many other things to occupy my time so don't go out of my way enough to find excuses to use it.

Optimising for fewest key strokes only makes sense transmitting to Pluto or beyond
  • Comment on Re^2: Quick script to check data logger data

Replies are listed 'Best First'.
Re^3: Quick script to check data logger data
by etj (Priest) on Aug 10, 2024 at 02:44 UTC
    I'd be interested to see a bit more concrete data and the problem you're solving, since I expect it would be easy for a journeyman like me to turn that into a couple of lines of PDL; that would then be a useful example I could then steal into some docs. Here, that would include a screenshot of how the Tk window actually looks, since PDL has some really good plotting stuff available.

      Turns out when I said "now we are starting to collect data", I spoke too soon. There are bugs in the software that drive the logger and, after months of back and forth between us and the manufacturer, we haven't been able to get the system to work correctly.

      The data consists of measurements in two axis from an induction magnetometer at 64Hz. The analysis looks for signals at fairly low frequency (7Hz maybe) that are related to some ionospheric phenomena I've forgotten about. Data files are a day's worth of data can come out to around 21MB. My test script looks at the first 16 seconds (1024 points) and plots a simple frequency/amplitude graph just to check that the data looks sane (basically that it has wiggly bits). The data below is fairly typical and shows a spike just below 8Hz - the interesting bit to the scientists. Absolute amplitude isn't very interesting for a quick assessment of the system.

      Optimising for fewest key strokes only makes sense transmitting to Pluto or beyond
        Now, the version that still uses Math::FFT, but PDL::Graphics::Simple only for visualisation (note the lack of needing to scale etc) - you need PDL, PDL::Graphics::Simple, and one underlying plotter, I recommend PDL::Graphics::Gnuplot:
        use strict; use warnings; use Math::FFT; use PDL; use PDL::Graphics::Simple; my $filePath = $ARGV[0] // "testData.txt"; my @data = do { open my $fh, $filePath or die "Can't open $filePath: $ +!"; local $/; split /\s+/, <$fh> }; my $fft = Math::FFT->new(\@data); my $spectrum = $fft->spctrm; shift @$spectrum; # Remove DC signal component line pdl($spectrum); print "ret> "; <STDIN>;
        Even if you prefer other methods of calculation, I hope you agree that the quick visualisations available with P:G:S are a game-changer. (The last line is just so the plot window doesn't get immediately closed because your program terminated - that's not a problem if you use a REPL like perldl)
        And finally (cue human-interest story), the version using PDL::FFTW3; unfortunately for me I had to learn about power spectra and then convert Math::FFT's spctrm code into a small PDL function:
        use strict; use warnings; use PDL; use PDL::FFTW3; use PDL::Graphics::Simple; my $filePath = $ARGV[0] // "testData.txt"; my @data = do { open my $fh, $filePath or die "Can't open $filePath: $ +!"; local $/; split /\s+/, <$fh> }; sub spectrum {my ($d)=@_; my $hn=$d->dim(0)-1; my $a2=$d->abs2/(($hn*2)**2); $a2->slice([1,$hn-1])*=2; $a2; } my $spectrum = spectrum(rNfft1 pdl(\@data)); $spectrum = $spectrum->slice('1:-1'); # Remove DC signal component line pdl($spectrum); print "ret> "; <STDIN>;
        If you're into IRC, then next time you feel you might benefit from using PDL but don't feel 100% clear on where to start, do feel free to join #pdl on irc.perl.org and ask questions! A lightweight way to do so is on the MetaCPAN page for various PDL modules, with the "Chat with us" link.
        What? Software/hardware with bugs??

        Once I'd emotionally recovered from dealing with that idea, I simplified your code into this (with your data above in a file you can give as an argument, or in the obvious default filename):

        use strict; use warnings; use Tk; use Tk::Canvas; use Math::FFT; my $filePath = $ARGV[0] // "testData.txt"; my @data = do { open my $fh, $filePath or die "Can't open $filePath: $ +!"; local $/; split /\s+/, <$fh> }; my $spectrum; my $fft = Math::FFT->new(\@data); $spectrum = $fft->spctrm; shift @$spectrum; # Remove DC signal component my $mw = MainWindow->new (-title => "Magnetometer Plotter"); my $canvas = $mw->Canvas (-height => 700, -width => 1024)->pack (); my $i = 0; $canvas->createLine( (map +(2 + 2 * $i++, $_), NormData($spectrum, 680)), -fill => 'blue' ); $mw->MainLoop; sub NormData { my ($data, $span) = @_; my ($min, $max); for my $datum (@$data) { $min //= $datum; $max //= $datum; $min = $datum if $min > $datum; $max = $datum if $max < $datum; } my $scale = $span / ($max - $min); map $span - ($_ - $min) * $scale + 10, @$data; }

        There are bugs in the software that drive the logger and, after months of back and forth between us and the manufacturer, we haven't been able to get the system to work correctly.

        You wouldn't happen to have more details on that sensor setup? Maybe a datasheet, interface description, or anything of that sort? Maybe i could help to come with an alternative solution of getting the data out of the sensor.

        PerlMonks XP is useless? Not anymore: XPD - Do more with your PerlMonks XP
        Also check out my sisters artwork and my weekly webcomics