Why is there noise for higher bandwidth intensity distributions?

Question

I have been trying to computationally evaluate the following E-field computationally:

It is basically the sum of waves of a given wavelength, weighted by the square root of a gaussian distribution for the wavelengths.

I compute it via Python by performing gauss quadrature integrals for each value of $x$, through the package scipy.integrate.quad function. The code is stated below:

# Imports
import numpy as np
import scipy as sp
from scipy import integrate

# Parameters
mu = 0.635 # mean wavelength
sigma = 0.01 # std dev of wavelength distribution
# wl is wavelength
x_ara = np.arange(0, 1.4, 0.01)

# Limits of Integration
lower, upper = mu - 4*sigma, mu+4*sigma
if lower < 0 :
    print('lower limit met')
    lower = 1e-15 # cannot evaluate sigma = 0 due to singularity for gaussian function

# Functions
def Iprofile_func(wl, mu, sigma):
    profile = np.exp(-( ((wl-mu) / (np.sqrt(2)*sigma))**2))
    return profile

def E_func(x_ara, wl, mu, sigma):
    return np.sqrt(Iprofile_func(wl, mu, sigma)) * np.cos(2*np.pi/wl * (x_ara))

# Computation
field_ara = np.array([])
for x in x_ara:
    def E(wl):
        return E_func(x, wl, mu, sigma)

    field = sp.integrate.quad(E, lower, upper)[0]
    field_ara = np.append(field_ara, field)

I fixed the value of $\mu$ = 0.635, and performed the same computation for two values of $\sigma$, $\sigma$ = 0.01 and $\sigma$ = 0.2. The arrays that I get, I have plotted below, the plot above is the wavelength distribution, while the plot below is the computed field array:

Why does noise appear in the computed field when the value of sigma increases?