2

我正在使用一个名为 emcee 的 Python 模块来对分发进行采样。我需要将一个 (37,100)(我分别命名为 Ntrig 和 Nsamp)数组events传递给下面的函数。

def mp(SNR2, *events):
    events = np.asarray(events).reshape((Ntrig,Nsamp))
    bessel = special.iv(0,np.sqrt(x*SNR2(event)))
    exp = np.exp(-0.5*(x+SNR2(event)))
    I = integrate.quad(lambda x: exp*bessel,0,SNRth**2)[0]
    return np.asarray([np.array[I for event in events[i]] for i in range(len(events))]).reshape(events.shape)

我不断收到错误:

ValueError: total size of new array must be unchanged

据我了解,*events会将events数组分解为 37*100 个单独的参数。不应该在我重塑数组的下一行将其放回 37 x 100 数组吗?

PS 在你问我为什么还要费心分解events成单独的参数之前 PS——模块需要这个才能工作,它不能接受一个数组。

完整回溯错误:

ValueError                                Traceback (most recent call last)
<ipython-input-17-c8e815326a69> in <module>()
----> 1 mp(SNR2,events)

<ipython-input-16-9f73f234c628> in mp(SNR2, *events)
      5 def mp(SNR2, *events):
      6     events = np.asarray(events).reshape((Ntrig,Nsamp))
----> 7     return np.asarray([np.array([integrate.quad(lambda x: np.exp(-0.5*(x+SNR2(event)))*special.iv(0,np.sqrt(x*SNR2(event))),0,SNRth**2)[0] for event in events[i]]) for i in range(len(events))]).reshape(events.shape)
      8 #    return integrate.quad(lambda x: 0.5*np.exp(-0.5*(x+SNR2(event)))*special.iv(0,np.sqrt(x*SNR2(event))),0,SNRth**2)[0]
      9 def pp(SNR2, *events):

/opt/local/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/scipy/integrate/quadpack.pyc in quad(func, a, b, args, full_output, epsabs, epsrel, limit, points, weight, wvar, wopts, maxp1, limlst)
    279         args = (args,)
    280     if (weight is None):
--> 281         retval = _quad(func,a,b,args,full_output,epsabs,epsrel,limit,points)
    282     else:
    283         retval = _quad_weight(func,a,b,args,full_output,epsabs,epsrel,limlst,limit,maxp1,weight,wvar,wopts)

/opt/local/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/scipy/integrate/quadpack.pyc in _quad(func, a, b, args, full_output, epsabs, epsrel, limit, points)
    343     if points is None:
    344         if infbounds == 0:
--> 345             return _quadpack._qagse(func,a,b,args,full_output,epsabs,epsrel,limit)
    346         else:
    347             return _quadpack._qagie(func,bound,infbounds,args,full_output,epsabs,epsrel,limit)

<ipython-input-16-9f73f234c628> in <lambda>(x)
      5 def mp(SNR2, *events):
      6     events = np.asarray(events).reshape((Ntrig,Nsamp))
----> 7     return np.asarray([np.array([integrate.quad(lambda x: np.exp(-0.5*(x+SNR2(event)))*special.iv(0,np.sqrt(x*SNR2(event))),0,SNRth**2)[0] for event in events[i]]) for i in range(len(events))]).reshape(events.shape)
      8 #    return integrate.quad(lambda x: 0.5*np.exp(-0.5*(x+SNR2(event)))*special.iv(0,np.sqrt(x*SNR2(event))),0,SNRth**2)[0]
      9 def pp(SNR2, *events):

<ipython-input-16-9f73f234c628> in SNR2(*events)
      1 def SNR2(*events):
----> 2     events = np.asarray(events).reshape((Ntrig,Nsamp))
      3     C = 5*np.pi**(-1.33333)*events**(1.66667)/(96*d**2)
      4     return C*integrate.quad(lambda f: f**(-2.3333)/S(f), 20, 1500, limit=1000)[0]
      5 def mp(SNR2, *events):

ValueError: total size of new array must be unchanged
4

1 回答 1

3

据我了解,events会将事件数组分解为 37*100 个单独的参数。

这不是真的。如果你打电话mp使用

mp(SNR2, events)

然后在里面mpevents将是一个 1 元素元组(arr,),其中arr是 (37, 100) 形数组。

如果你打电话mp使用

mp(SNR2, *events)

然后在里面mpevents将是一个 37 元素的元组,其中 37 个元素是 (37, 100) 形数组的 37 行。

如果你打电话mp使用

mp(SNR2, *events.flat)

然后在里面mpevents将是一个 37*100 元素的元组。


请注意,回溯的最后一节说:

<ipython-input-16-9f73f234c628> in SNR2(*events)
      1 def SNR2(*events):
----> 2     events = np.asarray(events).reshape((Ntrig,Nsamp))
      3     C = 5*np.pi**(-1.33333)*events**(1.66667)/(96*d**2)
      4     return C*integrate.quad(lambda f: f**(-2.3333)/S(f), 20, 1500, limit=1000)[0]
      5 def mp(SNR2, *events):

ValueError: total size of new array must be unchanged

因此,当 Python 在SNR2函数中时会引发错误。

由于在usingSNR2中被调用,并且是一个 (37,100) 形数组,因此变量 in是一个包含原始数组的 1 元素元组。那不是你想要的。mpSNR2(event)eventeventSNR2

修复代码的最简单方法是定义

def SNR2(events):  
    # no longer needed
    # events = np.asarray(events).reshape((Ntrig,Nsamp))

就像人们期望的那样传递事件。

但是,如果您不能更改 的签名SNR2,则必须使用

SNR2(*event.flat)

函数内部mp


参考:这里对解包操作符进行了很好的解释*,以及在定义函数和调用函数时如何使用语法。

于 2014-07-23T18:52:19.343 回答