8

I have a time-dependent signal.

I wish to plot its integration over time with time being x-axis and integration value being y-axis.

Is there any Python way of doing this?

To be more specific:

I have a time array, time, and a signal array, signal. They are of same dimension.

I need to integrate signal over time with scipy.integrate.trapz().

Instead of getting the final integral, I wish to see the integral varying as time passes.

lmjohns3
  • 7,422
  • 5
  • 36
  • 56
Sibbs Gambling
  • 19,274
  • 42
  • 103
  • 174

2 Answers2

10

Try using scipy.integrate.cumtrapz() instead :

plt.plot(time[:-1], scipy.integrate.cumtrapz(signal, x=time))
plt.show()

It computes an array containing the cumulative integral values.

http://docs.scipy.org/doc/scipy-0.10.1/reference/generated/scipy.integrate.trapz.html

lmjohns3
  • 7,422
  • 5
  • 36
  • 56
  • The x-axis is incorrect. Please note I have no fixed time interval. What I have is a time array. – Sibbs Gambling Aug 15 '13 at 06:05
  • Ah, in that case you can use your time array to indicate where to plot the y values. I've updated the answer. – lmjohns3 Aug 15 '13 at 06:28
  • exactly, when you run the code, here comes the problem. The x and y in the plot() are of different dimensions. Size of y is 20123, and x is 20124. How may I resolve that? – Sibbs Gambling Aug 15 '13 at 06:31
  • 1
    You can just trim off the last time value. This equates to aligning the integral values with the left side of the corresponding trapezoid. (You can trim off the first time value to align with the right side.) – lmjohns3 Aug 15 '13 at 06:32
  • How please? really new to Python – Sibbs Gambling Aug 15 '13 at 06:34
  • I already updated the answer ! `time[:-1]` will return all but the last element of the `time` array ; `time[1:]` will return all but the first element. This is called [slicing](http://stackoverflow.com/questions/509211/the-python-slice-notation). – lmjohns3 Aug 15 '13 at 06:37
5

A slightly better answer uses the optional "initial" argument. Here is a complete example:

import scipy.integrate as it
import numpy as np 
import matplotlib.pyplot as plt
t=np.linspace(0,1, 100)
y=t**2
y_int = it.cumtrapz(  y  ,  t, initial=0.0)  # y_int is same size as t
plt.plot(t, y_int)
plt.show()

This avoids the strange indexing like time[:-1]

wander95
  • 1,298
  • 1
  • 15
  • 22