使用多个值循环时得到不同的结果
Get different results when loop using more than one value
我正在尝试为运输编写一个光变曲线模拟,但是,当使用多个倾角值时(下面代码中的'ibound'),会出现一些奇怪的水平线最后的数字,但是如果每次只使用一个值就不会发生,e.g.ibound=[80,82] 是不行的,但是 ibound=[80] 或 ibound=[82] 给出了正确的结果。
import numpy as np
import math
import matplotlib.pyplot as plt
r_p=1.314*(7.149e7)
r_s=0.826*(6.957e8)
Area_p=np.pi*r_p**2.
P=1.3*86400.
omega=2.*np.pi/P
a=0.0226*(1.496e11)
def r_cen(t):
return a*np.sqrt( (np.sin(omega*t)**2) + (
(np.cos(math.radians(i))**2) *
(np.cos(omega*t)**2) ) )
def A_case2(rcen):
Aobs = 0.
nstep = 200
for step in range(0,nstep):
r=rcen-r_p+(float(step)*2.*r_p/float(nstep))
if r>r_s:
break
else:
theta = ((r**2.)+(rcen**2.)-(r_p**2.)) / (2.*r*rcen)
if theta > 1.:
theta = 1.
Aobs=Aobs+(2.*r*np.arccos(theta)*(r_p*2./float(nstep)))
return Aobs
LC=[]
phase = []
time = []
A=[]
ibound=[80,82]
for i in ibound:
for t in range(-int(P/5),int(P/5)):
rcen=r_cen(float(t))
phase.append(float(t)/P)
time.append(float(t))
if rcen>r_s+r_p:
A = 0.
elif r_s>rcen+r_p:
A = Area_p
else:
A = A_case2(rcen)
LC.append(1.-(A/(np.pi*r_s**2)))
yaxis = LC
xaxis = phase
plt.plot(xaxis,yaxis,'r',linewidth=1)
plt.xlim([-0.1,0.1])
plt.xlabel('phase')
plt.ylabel('normalised counts')
plt.show()
它是 plt.plot()
函数的神器。当您将所有结果放入同一个列表时,绘图仪将尝试使用一条连续线连接所有数据点。所以第一遍的结束数据点将连接到第二遍的第一个点。
如果更改绘图参数,您可以清楚地看到这一点,
plt.plot(xaxis, yaxis, 'or', linewidth=1)
这给出了,
如您所见,绘图仪没有尝试用连续线连接数据点,因此没有额外的水平线。
您的 xaxis
值在中途有一个不连续点,因此 xaxis
看起来像(假设步长为 0.1)[-0.2, -0.1...0.1, 0.2, -0.2, -0.1...0.1, 0.2]
。正如 DrBwts 所说,这就是 plt.plot
连接不连续点的位置,从而形成一条水平线。
如果您希望绘图仪仍然连接数据点(如果您希望图表在减小线宽后看起来不太起伏,这很有帮助),您可以替换代码中的最后一部分(来自 yaxis = LC
起)如下:
yaxis = LC
xaxis = phase
num_elements_per_list = len(phase)//2
x1 = xaxis[:num_elements_per_list]
x2 = xaxis[num_elements_per_list:]
y1 = yaxis[:num_elements_per_list]
y2 = yaxis[num_elements_per_list:]
plt.plot(x1, y1, 'r', linewidth=1)
plt.plot(x2, y2, 'r', linewidth=1)
plt.xlim([-0.1,0.1])
# This line is just to make the x axis a bit prettier
plt.xticks(np.arange(-0.1, 0.101, 0.05))
plt.xlabel('Phase')
plt.ylabel('Normalised Counts')
plt.show()
这会生成以下图:
python
ibound=[80,82,85] # the size of the list is not important
for i in ibound:
LC=[]
phase = []
time = []
A=[]
for t in range(-int(P/5),int(P/5)):
rcen=r_cen(float(t))
phase.append(float(t)/P)
time.append(float(t))
if rcen>r_s+r_p:
A = 0.
elif r_s>rcen+r_p:
A = Area_p
else:
A = A_case2(rcen)
LC.append(1.-(A/(np.pi*r_s**2)))
yaxis = LC
xaxis = phase
plt.plot(xaxis, yaxis, 'r', linewidth=0.5)
plt.xlim([-0.1, 0.1])
plt.xlabel('phase')
plt.ylabel('normalised counts')
plt.show()
如果你改变你的块,这将是你想要的。
为什么它结合了已经描述的所有行
Figure
我正在尝试为运输编写一个光变曲线模拟,但是,当使用多个倾角值时(下面代码中的'ibound'),会出现一些奇怪的水平线最后的数字,但是如果每次只使用一个值就不会发生,e.g.ibound=[80,82] 是不行的,但是 ibound=[80] 或 ibound=[82] 给出了正确的结果。
import numpy as np
import math
import matplotlib.pyplot as plt
r_p=1.314*(7.149e7)
r_s=0.826*(6.957e8)
Area_p=np.pi*r_p**2.
P=1.3*86400.
omega=2.*np.pi/P
a=0.0226*(1.496e11)
def r_cen(t):
return a*np.sqrt( (np.sin(omega*t)**2) + (
(np.cos(math.radians(i))**2) *
(np.cos(omega*t)**2) ) )
def A_case2(rcen):
Aobs = 0.
nstep = 200
for step in range(0,nstep):
r=rcen-r_p+(float(step)*2.*r_p/float(nstep))
if r>r_s:
break
else:
theta = ((r**2.)+(rcen**2.)-(r_p**2.)) / (2.*r*rcen)
if theta > 1.:
theta = 1.
Aobs=Aobs+(2.*r*np.arccos(theta)*(r_p*2./float(nstep)))
return Aobs
LC=[]
phase = []
time = []
A=[]
ibound=[80,82]
for i in ibound:
for t in range(-int(P/5),int(P/5)):
rcen=r_cen(float(t))
phase.append(float(t)/P)
time.append(float(t))
if rcen>r_s+r_p:
A = 0.
elif r_s>rcen+r_p:
A = Area_p
else:
A = A_case2(rcen)
LC.append(1.-(A/(np.pi*r_s**2)))
yaxis = LC
xaxis = phase
plt.plot(xaxis,yaxis,'r',linewidth=1)
plt.xlim([-0.1,0.1])
plt.xlabel('phase')
plt.ylabel('normalised counts')
plt.show()
它是 plt.plot()
函数的神器。当您将所有结果放入同一个列表时,绘图仪将尝试使用一条连续线连接所有数据点。所以第一遍的结束数据点将连接到第二遍的第一个点。
如果更改绘图参数,您可以清楚地看到这一点,
plt.plot(xaxis, yaxis, 'or', linewidth=1)
这给出了,
如您所见,绘图仪没有尝试用连续线连接数据点,因此没有额外的水平线。
您的 xaxis
值在中途有一个不连续点,因此 xaxis
看起来像(假设步长为 0.1)[-0.2, -0.1...0.1, 0.2, -0.2, -0.1...0.1, 0.2]
。正如 DrBwts 所说,这就是 plt.plot
连接不连续点的位置,从而形成一条水平线。
如果您希望绘图仪仍然连接数据点(如果您希望图表在减小线宽后看起来不太起伏,这很有帮助),您可以替换代码中的最后一部分(来自 yaxis = LC
起)如下:
yaxis = LC
xaxis = phase
num_elements_per_list = len(phase)//2
x1 = xaxis[:num_elements_per_list]
x2 = xaxis[num_elements_per_list:]
y1 = yaxis[:num_elements_per_list]
y2 = yaxis[num_elements_per_list:]
plt.plot(x1, y1, 'r', linewidth=1)
plt.plot(x2, y2, 'r', linewidth=1)
plt.xlim([-0.1,0.1])
# This line is just to make the x axis a bit prettier
plt.xticks(np.arange(-0.1, 0.101, 0.05))
plt.xlabel('Phase')
plt.ylabel('Normalised Counts')
plt.show()
这会生成以下图:
python
ibound=[80,82,85] # the size of the list is not important
for i in ibound:
LC=[]
phase = []
time = []
A=[]
for t in range(-int(P/5),int(P/5)):
rcen=r_cen(float(t))
phase.append(float(t)/P)
time.append(float(t))
if rcen>r_s+r_p:
A = 0.
elif r_s>rcen+r_p:
A = Area_p
else:
A = A_case2(rcen)
LC.append(1.-(A/(np.pi*r_s**2)))
yaxis = LC
xaxis = phase
plt.plot(xaxis, yaxis, 'r', linewidth=0.5)
plt.xlim([-0.1, 0.1])
plt.xlabel('phase')
plt.ylabel('normalised counts')
plt.show()
如果你改变你的块,这将是你想要的。 为什么它结合了已经描述的所有行
Figure