Files
very_short_lightning/2_get_flash_v2.py
2025-07-28 11:08:04 +08:00

92 lines
5.0 KiB
Python
Raw Permalink Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

import pandas as pd
from tqdm import tqdm
import os
df = pd.read_parquet("./data/flash/FLASH_2024.parquet")
'''
时间 微秒 纬度 经度 电流 回击 ... 误差椭圆长半轴方向 标志 陡度 时间差 RawID GUID
0 1/1/2024 03:16:40 5888618 25.850505 106.476598 -17.1 1 ... 118 0 0 0 1050961381,1050961380; 165612266
1 1/1/2024 03:42:33 7173002 25.684971 106.972500 73.1 1 ... 104 0 0 0 1050979230,1050979226,1050979231;1050979227,10... 165612269
2 1/1/2024 04:11:28 1344128 25.544373 106.873714 45.7 1 ... 87 0 0 0 1050991636,1050991653,1050991654;1050991651,10... 165612271
3 1/1/2024 04:37:00 2312000 26.191353 105.803525 -12.9 1 ... 118 0 0 0 1051001690,1051001680; 165612285
4 1/1/2024 05:22:22 1754035 26.291300 106.622184 -29.1 1 ... 0 0 0 0 1051017776,1051017786,1051017785; 165612299
'''
# 只保留需要的列并解析时间
df = df[['时间', '微秒', '纬度', '经度', '电流']].copy()
# 解析时间列,只保留日期和时分秒(不需要微秒)
# 先将'时间'列转换为datetime对象格式为日/月/年 时:分:秒
df['时间'] = pd.to_datetime(df['时间'], format='%d/%m/%Y %H:%M:%S')
df_station_id = pd.read_csv("station_id.csv")
EFID = df_station_id["EFID"].to_list()
Longitude = df_station_id["Longitude"].to_list()
Latitude = df_station_id["Latitude"].to_list()
for id,lat,lon in zip(EFID,Latitude,Longitude):
# 2. 定义目标点和范围 (与之前相同)
target_lon = lon
target_lat = lat
degree_range = 0.15
# 3. 筛选在指定经纬度范围内的点 (与之前相同)
df_filtered = df[
(df['经度'] >= target_lon - degree_range) &
(df['经度'] <= target_lon + degree_range) &
(df['纬度'] >= target_lat - degree_range) &
(df['纬度'] <= target_lat + degree_range)
].copy()
# 4. 核心逻辑修正:使用 resample 和 rolling
if df_filtered.empty:
# 输出未找到雷电数据的提示,并显示经纬度范围
lon_min: float = target_lon - degree_range # 计算经度下限
lon_max: float = target_lon + degree_range # 计算经度上限
lat_min: float = target_lat - degree_range # 计算纬度下限
lat_max: float = target_lat + degree_range # 计算纬度上限
print(f"{id}在指定范围内(经度: {lon_min:.5f}~{lon_max:.5f}, 纬度: {lat_min:.5f}~{lat_max:.5f})没有找到雷电数据。")
else:
# 将时间设置为索引
df_filtered.set_index('时间', inplace=True)
# 创建一个表示雷电事件的 Series每个事件计为 1
# 然后按5秒的间隔重采样并对每个间隔内的事件求和
# 这会得到一个时间序列其中索引是每5秒值是该5秒内发生的雷电次数
lightning_counts_per_5s = df_filtered.resample('5S').size()
# 使用滚动窗口计算未来30分钟的雷电数
# window='30min' 定义了窗口大小
# 为了实现“未来”30分钟的计数我们使用一个技巧
# 1. 将时间序列反转 (iloc[::-1])
# 2. 使用一个标准的向后看的滚动窗口
# 3. 再将结果反转回来
# 这样就等效于一个向前看的滚动窗口
# closed='right' 确保窗口包含其起始点(在反转前)
future_counts = lightning_counts_per_5s.iloc[::-1].rolling(window='30min', closed='right').sum().iloc[::-1]
# 上一步的结果只包含有雷电的5秒间隔。现在我们需要创建一个完整的时间范围。
# 定义我们关心的完整时间范围的开始和结束
start_time = df_filtered.index.min().floor('D') # 从数据第一天的零点开始
end_time = df_filtered.index.max().ceil('D') # 到数据最后一天的零点结束
# 创建一个从开始到结束每隔5秒一个点的完整时间索引
full_time_index = pd.date_range(start=start_time, end=end_time, freq='5S')
# 将我们计算出的结果填充到这个完整的时间索引上
# reindex 会自动对齐索引,没有值的地方会填充为 NaN我们再用 fillna(0) 填充为 0
final_result = future_counts.reindex(full_time_index).fillna(0).astype(int)
# 重命名列并计算flash_count的分布
result_df = final_result.reset_index()
result_df.columns = ['time', 'flash_count']
# 计算flash_count列中每个数值的数量
flash_count_distribution = result_df['flash_count'].value_counts().sort_index()
print("\nflash_count列中每个数值的数量")
print(flash_count_distribution)
# 您也可以将结果保存到文件
os.makedirs(f'data/preprocessed_flash/{id}/',exist_ok=True)
result_df.to_parquet(f'data/preprocessed_flash/{id}/flash_{id}.parquet')