INTERMAGNET地磁数据时间序列图绘制异常问题求助

问题背景 #

为开展研究，从INTERMAGNET下载了2013年1月1日至12月31日的THL（格陵兰Qaanaaq）地磁数据，尝试绘制X分量时间序列图，但结果与官网生成的图差异明显：

INTERMAGNET官网生成的图 #

自制图 #

使用的绘图代码如下：

import pandas as pd
import matplotlib.pyplot as plt

# 读取csv文件
df = pd.read_csv('data/thl_data.csv')

# 将日、月、年、时、分列转换为datetime
df['date'] = pd.to_datetime(df[['year', 'month', 'day', 'hour', 'minute']])

# 设置date为索引
df = df.set_index('date')

# 筛选2013年1月1日至12月31日的数据
start_date = '2013-01-01'
end_date = '2013-12-31'
df = df.loc[start_date:end_date]

# 绘制x值
plt.figure(figsize=(18,12))
plt.plot(df['x'], label='x')
plt.xlabel('Time')
plt.ylabel('X (nT)')
plt.title('THL (Qaanaaq (Thule), Greenland)')
plt.legend()
plt.show()

尝试过按日重采样取平均值，但绘制结果仍与未采样版本相似，无法匹配官网图。数据样例：

day,month,year,hour,minute,x,y,z
1,1,2013,0,1,26104,-31575,562205
1,1,2013,0,2,26105,-31584,562201
1,1,2013,0,3,26109,-31593,562197
1,1,2013,0,4,26115,-31597,562197
1,1,2013,0,5,26113,-31611,562190
1,1,2013,0,6,26112,-31605,562195
1,1,2013,0,7,26106,-31604,562195
1,1,2013,0,8,26111,-31594,562202
1,1,2013,0,9,26111,-31596,562203
1,1,2013,0,10,26116,-31598,562202
1,1,2013,0,11,26113,-31596,562203
1,1,2013,0,12,26114,-31599,562202
1,1,2013,0,13,26110,-31604,562201
1,1,2013,0,14,26114,-31598,562206
1,1,2013,0,15,26124,-31587,562211
...

解决思路与方法 #

1. 匹配INTERMAGNET的数据处理逻辑 #

INTERMAGNET官网展示的地磁图通常是去除长期趋势的扰动值，而非原始绝对值。原始数据包含年际/季节等缓慢变化的长期趋势，而官网图重点展示短周期地磁活动（比如地磁暴、日变化），因此需要：

• 计算X分量的长期基线（比如月均值或年均值），用原始X值减去基线得到扰动值ΔX，再绘制ΔX的时间序列。
• 示例代码：

# 计算月均值作为基线
df['x_monthly_mean'] = df['x'].resample('M').transform('mean')
# 计算扰动值
df['x_disturbance'] = df['x'] - df['x_monthly_mean']

# 绘制扰动值
plt.figure(figsize=(18,12))
plt.plot(df['x_disturbance'], label='ΔX')
plt.xlabel('Time')
plt.ylabel('ΔX (nT)')
plt.title('THL (Qaanaaq (Thule), Greenland) - X Component Disturbance')
plt.legend()
plt.show()

2. 使用小时均值而非原始分钟数据 #

INTERMAGNET的官方图一般采用小时均值来平滑高频噪声，而不是原始分钟数据：

# 按小时重采样取均值
df_hourly = df.resample('H').mean()

# 绘制小时均值的X分量（或其扰动值）
plt.figure(figsize=(18,12))
plt.plot(df_hourly['x'], label='Hourly Mean X')
plt.xlabel('Time')
plt.ylabel('X (nT)')
plt.title('THL (Qaanaaq (Thule), Greenland) - Hourly Mean X')
plt.legend()
plt.show()

3. 检查数据的坐标系统 #

部分地磁数据可能存在分量转换（比如从地磁坐标转换为地理坐标），确认下载的数据分量是否与官网展示的一致。INTERMAGNET通常提供的是地磁坐标下的X（北向）、Y（东向）、Z（垂直向下）分量，若你的数据是其他坐标系，需要进行转换。

4. 验证时间索引的正确性 #

确保pd.to_datetime生成的时间索引准确，比如确认时间序列无缺口或重复（你的数据样例中分钟从1开始属于正常的分钟级采样，无需调整）。

内容的提问来源于stack exchange，提问作者Nikhil Raghavendra