[AstroPy] Convolution with astropy

Stefano Roberto Soleti roberto.soleti at gmail.com
Thu Nov 7 14:34:36 EST 2019


Dear astropy users,

I am having a problem understanding how the convolution of two models is done in astropy. I want to convolve my function, centered around a certain value x, with a Gaussian with mean 0 and sigma > 0. However, when I try to plot the result of the convolution, the function is centered around ~x/2. I think I am missing something fundamental here but I couldn’t find anything apart from the basic example on astropy docs. 
Here it is my code snippet:

import numpy as np
import math
import matplotlib.pyplot as plt

from astropy.modeling.models import custom_model
from astropy.modeling import models
from astropy.convolution import convolve_models

from scipy.constants import alpha


@custom_model
def spectrum(x):
    me = 0.511
    eMax = 104.97

    E = np.sqrt(x*x + me*me)
    result = (1./eMax)*(alpha/(2*math.pi))*(np.log(4*E*E/me/me)-2.)*((E*E+eMax*eMax)/eMax/(eMax-E))
    
    result[result < 0] = 0
    return result

x = np.linspace(0,120,1000)
s = spectrum()
g1 = models.Gaussian1D(1, 0, 2)
conv = convolve_models(g1, s)

plt.plot(x, conv(x), 'k-')
plt.plot(x, s(x))


Thank you very much for your help!


More information about the AstroPy mailing list