2

A question on numerical integration (Trapezoid) on live data in Python-

Background: In real time I'm measuring a moving object which has a mean speed of 100 metres/min, and I'm sampling every 100ms for 60 seconds - therefore at the end of the minute I'll have 600 samples. My method of measurement is not accurate and has 10 metres variance per measurement. A typical graph after a minute interval would be:

enter image description here

Question: I know roughly that the object travels 100 metres in a minute, however when integrate over this period (just use Trapezoid for now) a I get an answer which is 60x more than expected, what am I doing wrong? I suspect it's the 'width' or deltaT = 100ms which is incorrect(?)

My python code is below - this is idiomatic i.e. not pythonic, for a reason; to mimic the real time measurements: 

# Trapezoidal rule integral
testData = []       # store all vel. measurements
width    = 100e-3   # sample every 100ms
total    = 0
currVel  = 0
prevVel  = 0
t = 0

while( t < 60 ):
    # take a live sample and save it
    currVel = np.random.normal(100,10,1)  
    testData.append( currVel )
    # only complete the integral after the second sample
    if( t >= 100e-3 ):
        total += width*(currVel+prevVel)/2
    # update previous flow and increment 
    prevVel = currVel
    t += 100e-3
#total /= 60  # fudge factor, why this factor out?
print "Total distance covered over 1 minute", total
Harry Lime
  • 2,167
  • 8
  • 31
  • 53
  • 1
    There's no reason to use trapezoidal here. If you think about it, you're adding half of a new sample at one time, then next cycle when that sample becomes prevVel, you're adding the other half. The final value in total is just the sum of all the samples, except for the endpoints. – DarenW Feb 23 '14 at 01:11
  • Thanks DarenW for the response but that doesn't answer my question: Why is it out by factor of 60, unless I've missed something? – Harry Lime Feb 23 '14 at 01:30
  • 1
    You have a units problem. Your sample distribution is m/min, try changing that to 1.67 m/s. Your loop is measured over 60 s, not 60 min. And what are you doing? Tracking people with radar? – Adam Burry Feb 23 '14 at 01:31

1 Answers1

1

An integral can be thought of as an area calculation. You essentially have:

(100 m/min) * (60 s)

and you are getting an answer of 6000 because the program has no representation of units. (The answer is 6000 meter-seconds per minute.) If you write your calculation this way instead, the mistake might be more obvious:

(100 m / 60 s) * (60 s)

Now the seconds will cancel out properly. See this thread for a discussion of Python unit libraries.

Community
  • 1
  • 1
Adam Burry
  • 1,904
  • 13
  • 20
  • Thanks Adam Burry, but I can't see where /60*60(?). That said I can see why you have to /60 (t=100ms/60); to accommodate the source rate of metres/min. I think this is sample rate conversion - http://en.wikipedia.org/wiki/Sample_rate_conversion – Harry Lime Feb 23 '14 at 13:04
  • And no, this is not related to sample rate conversion. – Adam Burry Feb 24 '14 at 02:16
  • OK Adam, so it must be a rate conversion then(?) forgive my naivety just wan to be more clear than 'Y/N type' responses. – Harry Lime Feb 28 '14 at 09:51