count.rb

#!/usr/bin/env ruby
require 'plw'

### Program
if ARGV[0].nil? or not File.exists? ARGV[0]
  $stderr.puts "Usage: #{$0} filename.plw\n  You must specify a PLW file to parse"
  Process.exit
end

$stdout.sync = true

filename = ARGV[0]

$stdout.puts "----------------------------------------"
$stdout.puts "|           PLW Peak Counter           |"
$stdout.puts "----------------------------------------"
$stdout.puts "Processing file: #{filename}\n\n"

# Load the PLW file
plw = PLW.new(filename)

chanspec = {}

if File.exists?(filename+".stats")
  $stdout.puts "Standard deviations already calculated. Skipping..."
  open(filename+".stats") do |f|
    chanspec = Marshal.load(f.read)
  end
else
  # Iterate through the first n points to determine the noise level (the average of the first n points, and the std.dev.)
  $stdout.puts "Calculating the standard deviation of the channel voltages\n"
  
  chanspec = [[0,0]] * plw.no_params
  n = 0
  last_print = nil
  while (sample = plw.getSample) != false
    i = -1
    chanspec.collect! { |chan| # Each plw channel has its own accumulator for the sum(data squared) and sum(data) squared needed for a std.dev.
      i += 1 # Channel Iterator
      [chan[0] + sample[:data][i], # Sum
      chan[1] + (sample[:data][i])**2] # Sum(squared)
    }
    pc = (n*100/plw.no_samples).round
    if pc != last_print
      ps = (pc * 0.7).round
      text = ">|#{pc}%"
      $stdout.write "\r|"+("="*[ps,70-text.length].min)+text+("-"*[(70 - ps - text.length),0].max)+"|"
      last_print = pc
    end
    n += 1
  end
  $stdout.write "\r - Complete, details saved                                              \n"
  chanspec.collect! { |chan|
    {
      :average => chan[0]/plw.no_samples,
      :std_dev => (chan[1]/plw.no_samples - (chan[0]/plw.no_samples)**2)**0.5
    }
  }
  open(filename+".stats","w") do |f|
    f.write Marshal.dump(chanspec)
  end
end

# Rewind to the first data point
plw.current_sample = 0

# We need to keep track of whether or not the previous point was marked as being a peak
# so that two consecutive points above the required level are still considered one peak
lastPeaks = [false] * plw.no_params
countPeaks = [0] * plw.no_params

$stdout.puts "Looking for peaks (voltage spikes outside 2 standard deviations of the mean)"

n = 0
last_print = nil
# Begin scanning for peaks
while (sample = plw.getSample) != false
  i = -1 # Channel Iterator
  # Measure to see if we're outsde the noise boundary, create a true/false array for the channels and also
  # add to the peak count if required
  lastPeaks = lastPeaks.collect{ |last|
    i += 1
    #
    is_peak = (sample[:data][i] - chanspec[i][:average]).abs > chanspec[i][:std_dev]*2
    if is_peak and not last
      countPeaks[i] += 1
    end
    is_peak
  }
  pc = (n*100/plw.no_samples).round
  if pc != last_print
    ps = (pc * 0.7).round
    text = ">|#{pc}%"
    $stdout.write "\r|"+("="*[ps,70-text.length].min)+text+("-"*[(70 - ps - text.length),0].max)+"|"
    last_print = pc
  end
  n += 1
end

$stdout.write "\r|============================ All Done! ===============================|\n"
p countPeaks