1
我試圖繪製數據以檢查我的代碼,我將所得繪圖與已經用Matlab生成的繪圖進行比較。然而,我遇到了幾個問題:在Python中繪製RINEX
- 一般來說,RINEX文件的解析工作,並且數據表示的一般模式看起來類似於Matlab腳本繪製的。但是,在數據放大時,數據會有小的偏差,例如在使用較小的時間序列時,例如繪製特殊的2小時而非24小時。在Matlab中,可以看到這個小的差異,並應用多項式擬合。但是對於Python圖(下面顯示的第一張圖),這兩小時期間的曲線看起來「平滑」,並且完全沒有偏離,就像在Matlab腳本中看到的那樣(第二個圖顯示藍色線作爲數據,相對於polyfit的紅線,因此,藍線在x = 9.4時表現出輕微的差異)。假定Matlab腳本是正確的,因爲這種偏差是由於地震活動暫時破壞了電離層。請參閱下面的情節:
第三個情節是在Matlab,其中這簡直是polyfit減去實時數據。
因此,不清楚該數據是如何繪製在Python腳本的座標軸上的,因爲數據看起來很平滑?也不如果我的代碼是錯誤的(見下文)並以某種方式「平滑」出的數據以某種方式:
#Calculating by looping through
for sv in range(32):
sat = self.obs_data_chunks_dataframe[sv, :]
#print "sat.index_{0}: {1}".format(sv+1, sat.index)
phi1 = sat['L1'] * LAMBDA_1 #Change units of L1 to meters
phi2 = sat['L2'] * LAMBDA_2 #Change units of L2 to meters
pr1 = sat['P1']
pr2 = sat['P2']
#CALCULATION: teqc Calculation
iono_teqc = COEFF * (pr2 - pr1)/1000000 #divide to make values smaller (tbc)
print "iono_teqc_{0}: {1}".format(sv+1, iono_teqc)
#PLOTTING
#Plotting of the data
plt.plot(sat.index, iono_teqc, label=‘teqc’)
plt.xlabel('Time (UTC)')
plt.ylabel('Ionosphere Delay (meters)')
plt.title("Ionosphere Delay on {0} for Satellite {1}.".format(self.date, sv+1))
plt.legend()
ax = plt.gca()
ax.ticklabel_format(useOffset=False)
plt.grid()
if sys.platform.startswith('win'):
plt.savefig(winpath + '\Figure_SV{0}'.format(sv+1))
elif sys.platform.startswith('darwin'):
plt.savefig(macpath + 'Figure_SV{0}'.format(sv+1))
plt.close()
從點1繼,多項式擬合代碼下面不按照我想要的方式運行,所以我在這裏忽略了一些東西。我認爲這與在x,y軸上使用的數據有關,但不能精確地確定是什麼。有人會知道我哪裏錯了嗎?
#Zoomed in plots if sv == 19: #Plotting of the data plt.plot(sat.index, iono_teqc, label=‘teqc’) #sat.index to plot for time in UTC plt.xlim(8, 10) plt.xlabel('Time (UTC)') plt.ylabel('Ionosphere Delay (meters)') plt.title("Ionosphere Delay on {0} for Satellite {1}.".format(self.date, sv+1)) plt.legend() ax = plt.gca() ax.ticklabel_format(useOffset=False) plt.grid() #Polynomial fitting coefficients = np.polyfit(sat.index, iono_teqc, 2) plt.plot(coefficients) if sys.platform.startswith('win'): #os.path.join(winpath, 'Figure_SV{0}'.format(sv+1)) plt.savefig(winpath + '\Zoom_SV{0}'.format(sv+1)) elif sys.platform.startswith('darwin'): plt.savefig(macpath + 'Zoom_SV{0}'.format(sv+1)) plt.close()
我的RINEX文件包含32顆衛星。然而試圖生成的曲線對於所有32時,我得到:
IndexError:指數31超出範圍爲軸線0與尺寸31
變更下方至31個解決了這個部分中,只有代碼不包括第32顆衛星。我想也爲衛星32繪圖。解析和格式化數據的功能如下:
def read_obs(self, RINEXfile, n_sat, sat_map):
obs = np.empty((TOTAL_SATS, len(self.obs_types)), dtype=np.float64) * np.NaN
lli = np.zeros((TOTAL_SATS, len(self.obs_types)), dtype=np.uint8)
signal_strength = np.zeros((TOTAL_SATS, len(self.obs_types)), dtype=np.uint8)
for i in range(n_sat):
# Join together observations for a single satellite if split across lines.
obs_line = ''.join(padline(RINEXfile.readline()[:-1], 16) for _ in range((len(self.obs_types) + 4)/5))
#obs_line = ''.join(padline(RINEXfile.readline()[:-1], 16) for _ in range(2))
#while obs_line
for j in range(len(self.obs_types)):
obs_record = obs_line[16*j:16*(j+1)]
obs[sat_map[i], j] = floatornan(obs_record[0:14])
lli[sat_map[i], j] = digitorzero(obs_record[14:15])
signal_strength[sat_map[i], j] = digitorzero(obs_record[15:16])
return obs, lli, signal_strength
def read_data_chunk(self, RINEXfile, CHUNK_SIZE = 10000):
obss = np.empty((CHUNK_SIZE, TOTAL_SATS, len(self.obs_types)), dtype=np.float64) * np.NaN
llis = np.zeros((CHUNK_SIZE, TOTAL_SATS, len(self.obs_types)), dtype=np.uint8)
signal_strengths = np.zeros((CHUNK_SIZE, TOTAL_SATS, len(self.obs_types)), dtype=np.uint8)
epochs = np.zeros(CHUNK_SIZE, dtype='datetime64[us]')
flags = np.zeros(CHUNK_SIZE, dtype=np.uint8)
i = 0 #ggfrfg
while True:
hdr = self.read_epoch_header(RINEXfile)
if hdr is None:
break
epoch_time, flags[i], sats = hdr
#epochs[i] = np.datetime64(epoch_time)
epochs[i] = epoch_time
sat_map = np.ones(len(sats)) * -1
for n, sat in enumerate(sats):
if sat[0] == 'G':
sat_map[n] = int(sat[1:]) - 1
obss[i], llis[i], signal_strengths[i] = self.read_obs(RINEXfile, len(sats), sat_map)
i += 1
if i >= CHUNK_SIZE:
break
return obss[:i], llis[:i], signal_strengths[:i], epochs[:i], flags[:i]
def read_data(self, RINEXfile):
obs_data_chunks = []
while True:
obss, _, _, epochs, _ = self.read_data_chunk(RINEXfile)
epochs = epochs.astype(np.int64)
epochs = np.divide(epochs, float(3600.000))
if obss.shape[0] == 0:
break
obs_data_chunks.append(pd.Panel(
np.rollaxis(obss, 1, 0),
items=['G%02d' % d for d in range(1, 33)],
major_axis=epochs,
minor_axis=self.obs_types
).dropna(axis=0, how='all').dropna(axis=2, how='all'))
self.obs_data_chunks_dataframe = obs_data_chunks[0]
有什麼建議嗎?
乾杯,pymat。