黑体问题和卡方最小化误差

问题描述

我是 Python 的大一新生，我正在尝试编写一个代码来计算黑体分布问题的正确温度。 有问题的物体是宇宙，我正在使用来自 FIRAS CMB 单极谱 (you can get them here) 的数据。

问题是：当我用卡方最小化时，卡方的值每次总是返回相同的值。

我不知道问题出在哪里或出在哪里，此时，我厌倦了尝试不同的事情。

顺便说一句，我只想将卡方的每个值存储到一个数组中，这就是我所需要的。

import numpy as n
import scipy.constants as scy
import pandas as p 
a=p.read_csv('spec-1-mod.csv') 
b=a.loc[:,"column1"]
nn=a.loc[:,"column2"]
n2n=a.loc[:,"column4"]
def disPlanckTeoric(nu,Tin):
    jj=2*scy.h*scy.c*(nu**3)/scy.c**2
    jk=(scy.h*nu)/scy.k*Tin
    return jj*(1/(n.exp(jk)-1))
T=300
Tsize=n.arange(T,-1)
expe=n.array(nn)
otra=n.array(n2n) 
for j in range (len(Tsize)):
    teo=n.array(disPlanckTeoric(b,T)) 
    chi=n.sqrt(sum((((expe-teo))/otra)**2))
    chiM=n.empty((2,len(Tsize)))
    for k in range(len(Tsize)):
        chiM[:,k]=chi
    T=T-0.1
#print(chiM) #Just to see that the values are all the same

解决方法

您的代码中有一些输入错误。例如普朗克分布，分子中没有 c，在 jk=(scy.h*nu)/scy.k*Tin 中，您应该将 scy.k*Tin 放在括号中：jk=(scy.h*nu)/(scy.k*Tin)。在您的代码 chi=n.sqrt(sum((((expe-teo))/otra)**2)) 中，正确的部分是相同的数字，因为 sum 始终为 j 的每个值返回一个数字。

但是我更改了您的代码。我没有使用 pandas 和 numpy。我什至不确定我的代码绘制的曲线是你想要的。但请检查一下。

from math import *
from scipy.constants import *
import matplotlib.pyplot as plt

lines = open('firas_monopole_spec_v1.txt','r').readlines()
lines = [ line for line in lines if line[0] != '#' ] 
column0 = []
column1 = []
column2 = []
column3 = []
column4 = []
for line in lines:
    line_list = [ element for element in line.split(' ') if element != '' ]
    line_list[-1] = line_list[-1][:-1]
    column0.append( line_list[0] )
    column1.append( line_list[1] )
    column2.append( line_list[2] )
    column3.append( line_list[3] )
    column4.append( line_list[4] )

def disPlanckTeoric(nu,Tin):
    jj=2*h*(nu**3)/c**2
    jk=(h*nu)/(k*Tin)
    return jj*(1/(exp(jk)-1))

T=300
b = column0
expe = column1
otra = column3
teo = []
chi = []
tmp = []
T_arr = []
for j in range(0,len(b)):
    T_arr.append(j+1)
    for i in range(len(b)):
        teo.append( disPlanckTeoric(float(b[i]),T_arr[j]) )
        tmp.append( ( ( float(expe[i])-float(teo[i]) ) / float(otra[i]) )**2 )
    chi.append( sqrt( sum( tmp) ) )

plt.plot(chi,T_arr)
plt.show()