2013-05-09 2 views
0

Я написал код, который позволяет мне принимать данные из файлов DICOM и разделять данные на отдельные сигналы FID.«аргумент должен быть вызываемой функцией» при использовании scipy.integrate.quad

from pylab import * 
import dicom 
import numpy as np 
import scipy as sp 

plan=dicom.read_file("1.3.46.670589.11.38085.5.22.3.1.4792.2013050818105496124") 
all_points = array(plan.SpectroscopyData) 
cmplx_data = all_points[0::2] -1j*all_points[1::2] 
frames = int(plan.NumberOfFrames) 
fid_pts = len(cmplx_data)/frames 
del_t = plan.AcquisitionDuration/(frames * fid_pts) 

fid_list = [] 
for fidN in arange(frames): 
    offset = fidN * fid_pts 
    current_fid = cmplx_data[offset:offset+fid_pts] 
    fid_list.append(current_fid) 

Я бы сейчас хотел бы количественно оценить некоторые из этих данных так, после применения преобразования Фурье и сдвига Я пытался использовать функцию четверной SciPy и возникла следующая ошибка:

spec = fftshift(fft(fid_list[0]))) 
sp.integrate.quad(spec, 660.0, 700.0) 



error          Traceback (most recent call last) 
/home/dominicc/Experiments/In Vitro/Glu, Cr/Phantom 1: 10mM Cr, 5mM Glu/WIP_SV_PRESS_ME_128TEs_5ms_spacing_1828/<ipython-input-107-17cb50e45927> in <module>() 
----> 1 sp.integrate.quad(fid, 660.0, 700.0) 

/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.pyc in quad(func, a, b, args, full_output, epsabs, epsrel, limit, points, weight, wvar, wopts, maxp1, limlst) 
243  if type(args) != type(()): args = (args,) 
244  if (weight is None): 
--> 245   retval = _quad(func,a,b,args,full_output,epsabs,epsrel,limit,points) 
246  else: 
247   retval = _quad_weight(func,a,b,args,full_output,epsabs,epsrel,limlst,limit,maxp1,weight,wvar,wopts) 

/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.pyc in _quad(func, a, b, args, full_output, epsabs, epsrel, limit, points) 
307  if points is None: 
308   if infbounds == 0: 
--> 309    return _quadpack._qagse(func,a,b,args,full_output,epsabs,epsrel,limit) 
310   else: 
311    return _quadpack._qagie(func,bound,infbounds,args,full_output,epsabs,epsrel,limit) 

error: First argument must be a callable function. 
МОГ

кто-нибудь, пожалуйста, предложите способ сделать этот объект доступным? После прочтения What is a "callable" in Python? я все еще не совсем понятен. После прочтения этого я попробовал

def fid(): 
    spec = fftshift(fft(fid_list[0])) 
    return spec 

Который возвратил ту же ошибку.

Любая помощь была бы принята с благодарностью, спасибо.

+0

Разве это не функция, которая принимает аргумент? Как еще это можно интегрировать? –

+0

Этот вопрос не имеет смысла. Вы не можете использовать 'scipy.integrate.quad' для интеграции массива. Он объединяет * функцию *, которая принимает аргумент (например, 'sin (x)'). Используйте что-то вроде 'scipy.integrate.simps' для интеграции массива значений функций. – talonmies

ответ

2

Поскольку вы интегрируете функцию, определенную только на сетке (т. Е. Массив чисел), вам нужно использовать одну из подпрограмм из «Integration, given fixed samples».

Фактически, первый аргумент quad()must be a function, то, что вы можете назвать («вызываемым»). Например, вы можете сделать:

>>> from scipy.integrate import quad 
>>> def f(x): 
...  return x**2 
... 
>>> quad(f, 0, 1) 
(0.33333333333333337, 3.700743415417189e-15) 

Это может быть сделано вместо этого путем отбора проб f():

>>> from scipy.integrate import simps 
>>> x_grid = numpy.linspace(0, 1, 101) # 101 numbers between 0 and 1 
>>> simps(x_grid**2, dx=x_grid[1]-x_grid[0]) # x_grid**2 is an *array*. dx=0.01 between x_grid values 
0.33333333333333337 

Ваша ситуация, как во втором примере: интегрировать выборочную функцию, поэтому естественно используйте один или cumtrapz(), simps() или romb().

Смежные вопросы