import astropy.io.fits as ap
    import numpy as np
    import  matplotlib 
    import matplotlib.pyplot as plt
    from scipy.interpolate import interp1d
    from scipy.ndimage import gaussian_filter
    import math
    from scipy.special import sph_harm
    from scipy.special import legendre
    import scipy.integrate

    def legendre(x,n):
     leg = legendre(n)
     P_n = leg(x)
     return P_n




    hdu1 = ap.open("magnetogram.fits")
    hdu1.info()  # so only one hdu list is there in this magnetogram fits file with 360X180
    data = hdu1[0].data
#print(data[0:,359])  # the first run from 0 to 179 and other run from 0 to 359



    theta_old = np.arccos(np.linspace(-1,1,num=180))
    theta_new = np.linspace(180*(3.14/180),180)

    Br_new = np.empty((180,360),float) #empty array of the same dimesnion as Br magnetogram

    for i in range(0,360):
     f = interp1d(theta_old,data[0:,i],kind="cubic")
     Br_new[0:,i] = f(theta_new)

    Br_new = gaussian_filter(Br_new,sigma=[2,2])

    print(data.shape,theta_old.shape)
    print(Br_new.shape)
    print(Br_new)

    plt.figure()
#x =  theta_old
#y = np.linspace(0,360,361)
#Z = data
#plt.imshow(Z,interpolation= "bilinear")
#plt.show()

    x_new = (np.linspace(3.14,181))
    y = np.linspace(0,361)
    Z = Br_new
    plt.imshow(Z,interpolation= "bilinear")
    plt.show()
    print("nnn")

    CmlCoefficients = []

    for l in range(0,3):
     #a = []
     for m in range(-l,l+1):
        freal = lambda Theta,Phi: Br_new[int(Theta),int(Phi)] * math.sin((3.14 / 180) * Theta) * 
               (np.real(
            scipy.special.sph_harm(-m,l,(3.14 / 180) * Phi,(3.14 / 180) * Theta)))  # legendre polynomial handle

        fimaginary = lambda Theta,int(Phi)] * math.sin((3.14 / 180) * Theta) * 
               (np.imag(
            scipy.special.sph_harm(-m,(3.14 / 180) * Theta)))  # legendre 
         polynomial handle

        answerreal = (scipy.integrate.dblquad(freal,lambda Theta: 0,lambda Theta: 180))
        answerimag = (scipy.integrate.dblquad(fimaginary,lambda Theta: 180))
        CmlCoefficients.append((l,m,complex(answerreal[0],answerimag[0])))
        print(answerreal,answerimag)

        print(CmlCoefficients)

我正在尝试对一个函数F（theta，Phi）进行积分，其中被积函数依赖于具有某些特殊数学函数（球谐函数）的产品中形状为（180，360）的二维数据数组。现在，我的代码可以工作了，并评估了l和m的每个值的积分，但是积分中的误差太大了。例如，对于l = 0和m = 0，积分的值为-5113.672846048103，误差为257.68725625979823

有什么办法可以减少错误？