wordpress 给文章添加幻灯武汉网络优化知名乐云seo
文章目录
- 前言
- 一、基本说明
- 二、代码
前言
该 MATLAB 代码用于从 GeoTIFF 文件中提取基于特定地理位置(经纬度)和日期的某个点的相关数据。代码首先读取一个包含事件数据(日期、经纬度)的 Excel 文件,然后根据日期和位置尝试从存储滑坡风险数据的 GeoTIFF 文件中提取相应的数值。为了加快数据提取过程,代码使用了并行处理。
一、基本说明
数据加载:从 Excel 文件中读取包含日期、经纬度的事件数据。
路径设置:指定 GeoTIFF 文件存储的文件夹路径。
GeoTIFF 数据提取:针对每个事件的日期和位置:
构建文件名并检查文件是否存在。
若文件存在,计算指定经纬度对应的行列索引并提取相应数值。
若坐标超出边界,则赋值为 NaN。
并行处理:使用 parfor 来加速处理每个日期-位置对。
结果存储:将提取的数值存储在数组 gpmresult 或 mswepresult 中。
该代码使得从大规模数据集中高效提取滑坡风险数据成为可能,并且具备处理缺失或越界数据的鲁棒性。
二、代码
% 清空环境、命令窗口,并关闭所有打开的图形窗口
clear all;
clc;
close all;% 加载滑坡事件数据,包括日期、纬度和经度
data = readtable('E:\data.xlsx');
longitude = data{:, 3}; % 提取经度数据
latitude = data{:, 2}; % 提取纬度数据
dates = data{:, 1}; % 提取日期,格式为 'yyyy/MM/dd'% 设置 GeoTIFF 文件夹路径
tif_folder = 'E:\';% 初始化结果数组
gpmresult = [];% 使用并行处理遍历每个数据点
parfor i = 1:length(dates)% 检查日期是单元格格式还是字符数组if iscell(dates)date_str = dates{i};elsedate_str = dates(i, :);end% 将日期字符串转换为 datetime 对象date_obj = datetime(date_str, 'InputFormat', 'yyyy/MM/dd');year_folder = datestr(date_obj, 'yyyy'); % 从日期中提取年份作为文件夹名称% 获取当前点的经纬度lon = longitude(i);lat = latitude(i);% 使用 'yyyymmdd' 格式构建 GeoTIFF 文件路径tif_file = fullfile(tif_folder, year_folder, sprintf('%s.tif', datestr(date_obj, 'yyyymmdd')));% 检查 GeoTIFF 文件是否存在if isfile(tif_file)% 获取 GeoTIFF 文件的元数据并读取数据info = geotiffinfo(tif_file);[data, R] = readgeoraster(tif_file);% 计算指定纬度和经度的行列索引row = ceil((max(R.LatitudeLimits) - lat) / R.CellExtentInLatitude);col = ceil((lon - min(R.LongitudeLimits)) / R.CellExtentInLongitude);% 检查计算的行列是否越界if row <= 0 || col <= 0 || row > size(data, 1) || col > size(data, 2)value30 = NaN; % 如果越界则设置为 NaNelsevalue30 = data(row, col); % 如果在范围内则提取数据值end% 将提取的值添加到结果数组中gpmresult = [gpmresult; value30];else% 如果文件不存在则输出提示信息fprintf('File not found: %s\n', tif_file);end
end
disp('GPM well done');% 第
% 以下部分与第一个代码块类似,但用于另一个数据集
% 使用日期格式 'yyyydoy' (一年中的第几天)% % 设置 MSWEP GeoTIFF 文件夹路径
% tif_folder = 'H:\';
% mswepresult = [];
%
% % 使用并行处理遍历每个点
% parfor i = 1:length(dates)
% % 检查日期格式并转换
% if iscell(dates)
% date_str = dates{i};
% else
% date_str = dates(i, :);
% end
%
% date_obj = datetime(date_str, 'InputFormat', 'yyyy/MM/dd');
% year_folder = datestr(date_obj, 'yyyy');
% day_of_year = day(date_obj, 'dayofyear'); % 获取一年中的第几天
%
% lon = longitude(i);
% lat = latitude(i);
% tif_file = fullfile(tif_folder, year_folder, sprintf('_%d%03d.tif', year(date_obj), day_of_year));
%
% if isfile(tif_file)
% info = geotiffinfo(tif_file);
% [data, R] = readgeoraster(tif_file);
%
% row = ceil((max(R.LatitudeLimits) - lat) / R.CellExtentInLatitude);
% col = ceil((lon - min(R.LongitudeLimits)) / R.CellExtentInLongitude);
%
% if row <= 0 || col <= 0 || row > size(data, 1) || col > size(data, 2)
% value30 = NaN; % 如果越界则设置为 NaN
% else
% value30 = data(row, col); % 如果在范围内则提取数据值
% end
%
% mswepresult = [mswepresult; value30];
% else
% fprintf('File not found: %s\n', tif_file);
% end
% end
% disp('MSWEP well done');