我在读取一个nc数据的时候,想通过读取一系列用户自定义的经纬度坐标来获取计算出来的温度数据 。但是,我遇到一个问题:nc数据本身没有经纬度信息!数据本身的情况如下:
nc数据储存的温度数据TBOT信息有三维:time=124, lat=360, lon=720。因为储存了每天每6小时的温度数据,一天有4个6小时,所以一个月一共有31*4=124,即time=124。 lat代表纬度,即从南纬-90°到北纬90°,每0.5°为一个步长储存温度数据,所以lat=(90+90)/0.5=360。 lon代表经度,即从0°到360°,每0.5°为一个步长储存温度数据,所以lat=360/0.5=720。
因为三维数据打不开,我对数据进行了重塑压缩,形成259200x124 double的矩阵t_2d,并进行了后续一系列计算。259200=720*360,相当于用一系列经纬度索引确定了一个位置,这个大家可以理解吧。
问题来了!259200这每一行代表的是地球上每一点的位置的温度值,但是这个位置是从哪里开始,朝着什么方向,在哪里结束呢?我猜测是从(0.25 ,-89.75)这个地方开始到(359.75,89.75)结束。所以,我想通过计算输入的经纬度相对于起始位置(0.25 ,-89.75)的偏移量来找到指定经纬度值对应的索引,然后将这个索引的经纬度相乘,得到的最终数值就是该位置对应259200中的行数。可是我发现这似乎不对。我发现一些地方的温度与其季节明显不符,一看就知道错的。所以这个问题,应该怎么解决
clc
clear all;
% 设置循环参数
start_month = 7;
end_month = 8;
% 存储每日和每月各温度指标的变量
avg_t_daily = [];
max_temp_daily= [];
min_temp_daily= [];
avg_t_monthly = [];
max_t_monthly = [];
min_t_monthly = [];
temperature=[];
studyaera=[];
LA=[];
% 循环读取和处理多个nc图像
for month = start_month:end_month
% 构建nc图像文件名
file_name = sprintf('F:\\TPHWL6Hrly\\clmforc.cruncep.V7.c2016.0.5d.TPQWL.2015-%02d.nc', month);
% 读取温度数据
t = ncread(file_name, 'TBOT');
LONGXY= ncread(file_name, 'LONGXY');
LATIXY= ncread(file_name, 'LATIXY');
% 获取变量 t 的尺寸信息(依次为:经度、纬度、时间),并将其分别存储在自定义变量 num_lon、num_lat 和 num_time
[num_lon, num_lat, num_time] = size(t);
% 重塑矩阵为二维
t_2d = reshape(t, [], num_time);
%% 求每天的气温情况,包括平均气温,日高温,日低温,365个值
% 计算365天每天的平均气温
num_days = floor(num_time / 4); % 总共的天数,向下取整
avg_t_daily_month = zeros(size(t_2d, 1), num_days); % 存储每天的平均温度
for i = 1:num_days
start_index = (i - 1) * 4 + 1;
end_index = i * 4;
avg_t_daily_month(:, i) = mean(t_2d(:, start_index:end_index), 2) - 273.15;
end
% 将每日平均温度存储到整体变量中
avg_t_daily = [avg_t_daily avg_t_daily_month];
% 计算365天每天的最高温度
max_temp_daily_month = zeros(size(t_2d, 1), num_days); % 存储每天的最高温度
for i = 1:num_days
start_index = (i - 1) * 4 + 1;
end_index = i * 4;
max_temp_daily_month(:, i) = max(t_2d(:, start_index:end_index), [], 2) - 273.15;
end
% 将365天每天最高温度存储到整体变量中
max_temp_daily = [max_temp_daily max_temp_daily_month];
% 计算365天每天的最低温度
min_temp_daily_month = zeros(size(t_2d, 1), num_days); % 存储每天的最高温度
for i = 1:num_days
start_index = (i - 1) * 4 + 1;
end_index = i * 4;
min_temp_daily_month(:, i) = min(t_2d(:, start_index:end_index), [], 2) - 273.15;
end
% 将计算365天每天的最低温度存储到整体变量中
min_temp_daily = [min_temp_daily min_temp_daily_month];
end
%% 求每月的气温情况,包括平均气温,月高温,月低温,12个值
% 计算每月平均气温
avg_t_mons_month = mean(t_2d, 2) - 273.15; % 沿纬度方向求平均值
% 将每个月的平均温度存储到整体变量中
avg_t_monthly = [avg_t_monthly avg_t_mons_month];
% 计算每月最高温
max_t_mons_month = max(t_2d, [], 2) - 273.15; % 沿纬度方向求平均值
% 将每个月的最高温存储到整体变量中
max_t_monthly = [max_t_monthly max_t_mons_month];
% 计算每月最低温
min_t_mons_month = min(t_2d, [], 2) - 273.15; % 沿纬度方向求平均值
% 将每个月的最低温存储到整体变量中
min_t_monthly = [min_t_monthly min_t_mons_month];
%% 求自定义地区的气温情况
longitudes = [113.17, 47.49, 110.75];
latitudes = [23.8, 30.30, 22.75];
% 定义纬度和经度的范围
lat_range = [-89.75, 89.75];
lon_range = [0.25, 359.75];
% 计算纬度和经度的步长
lat_step = 0.5;
lon_step = 0.5;
% 找到指定经纬度值对应的索引
lat_indices = floor((latitudes - lat_range(1)) / lat_step) + 1;
lon_indices = floor((longitudes - lon_range(1)) / lon_step) + 1;
% 经纬度值组合最终索引
result = lat_indices + (lon_indices - 1) * (lat_range(2) - lat_range(1) + 1) / lat_step;
% 根据最终索引循环获取对应的 t_2d 值
t_2d_studyarea = t_2d(result, :)-273.15;
avg_t_daily_studyarea = avg_t_daily(result, :)
min_temp_daily_studyarea = min_temp_daily(result, :)
max_temp_daily_studyarea = max_temp_daily(result, :)
avg_t_mons_month_studyarea = avg_t_mons_month(result, :)
max_t_mons_month_studyarea = max_t_mons_month(result, :)
min_t_mons_month_studyarea = min_t_mons_month(result, :)