logoChat01.ai

下列代码的作用是什么 %% 基础参数设置 years = 1980:2024; ...

Created on: May 18, 2025

Answered using O4 Mini High by Chat01

Question

下列代码的作用是什么
%% 基础参数设置
years = 1980:2024; % 分析年份
time_bins = 0:6:24; % 地方时划分
seasons = {'mam', 'jja', 'son', 'djf', 'annual'};
season_info = struct('months', {[3,4,5], [6,7,8], [9,10,11], [12,1,2]});

%% 路径配置
data_root = '/workdir/tanjy/data/era5/cp/';
file_pattern = 'cp_%d%02d.nc';

%% 地理参数
load('/workdir/tanjy/data/landmask_era5.mat');
load('/workdir/tanjy/data/coastlines.mat');
lat_vector = 90:-0.25:-60;
lon_vector = 0:0.25:359.75;
lats = length(lat_vector);
lons = length(lon_vector);

%% 初始化数据存储（四维：纬度×经度×时段×年份）
pixel_cp_mm = cell(5, 1);
for s = 1:5
pixel_cp_mm{s} = zeros(lats, lons, 4, length(years));
end

%% 主处理流程
years_data = 1980:2025;
for y = 1:length(years_data)
current_year = years_data(y);
for month = 1:12
if (current_year == 1980 && month <3) || (current_year == 2025 && month > 2)
continue;
end

text
    % 读取nc文件
    ncfile = fullfile(data_root, sprintf(file_pattern, current_year, month));
    fprintf('Processing %d%02d (%d/%d)\n', current_year, month, y, length(years_data));

    % 读取时间信息
    time_sec = ncread(ncfile, 'valid_time');
    time_dt = datetime(time_sec, 'ConvertFrom', 'epochtime', 'Epoch', '1970-01-01');

    % 三维数据读取
    start = [1, 1, 1];                            % [经度起始, 纬度起始, 时间起始]
    count = [1440, lats, Inf];                    % [读取经度数, 读取纬度数, 时间数]
    precip = ncread(ncfile, 'cp', start, count); % 维度原始为[经度, 纬度, 时间]
    precip = permute(precip, [2,1,3]);            % 维度调整为[纬度, 经度, 时间]

    % 时间步循环
    for t = 1:size(precip, 3)
        dt = time_dt(t);
        utc_hour = hour(dt);

        % 计算地方时
        local_time = utc_hour + lon_vector(:, :)/15;
        local_time = mod(local_time, 24);
        if local_time == 0
            local_time = 24;
        end

        % 时段分类
        for h = 1:4
            time_mask = (local_time > time_bins(h)) & (local_time <= time_bins(h+1));
            precip_h = precip(:,:,t)*1000 .* time_mask .* landmask_era5(1:601, :); % 应用陆地掩膜

            % 排除nan值
            if nnz(isnan(precip_h(:,:,:))) ~= 0
                continue;
            end

            % 确定季节归属
            for s_check = 1:4
                if any(month == season_info(s_check).months)
                        
                    % 计算季节年
                    if s_check == 4
                        if month == 12
                            season_year = current_year;
                        else
                            season_year = current_year - 1;
                        end
                    else
                        season_year = current_year;
                    end
                    if season_year >= 1980 && season_year <= 2024
                        yr_idx = season_year - 1979;
                        pixel_cp_mm{s_check}(:, :, h, yr_idx) = ...
                        pixel_cp_mm{s_check}(:, :, h, yr_idx) + precip_h;
                    end
                    break;
                end
            end
        end
    end

    % 文件处理完成清理
    clear ncfile time_sec time_dt precip utc_hour local_time time_mask precip_h 
    drawnow; % 强制刷新内存
end
fprintf('└─ 年份%d处理完成\n\n', current_year);

end

%% 全年降水累加
pixel_cp_mm{5} = sum(cat(5, pixel_cp_mm{1:4}), 5);

%% 百分比转换
pixel_cp_pct = cell(5, 1);
for s = 1:5
total = sum(pixel_cp_mm{s}, 3);
total(total == 0) = eps; % 防除零处理
pixel_cp_pct{s} = bsxfun(@rdivide, pixel_cp_mm{s}, total);
end

%% 计算斜率
pixel_kcp_mm = cell(5, 1);
pixel_kcp_pct = cell(5, 1);

% 准备回归参数
X = [ones(length(years),1), years(:)-mean(years)]; % 中心化设计矩阵

for s = 1:5

text
% 数据维度重排
Y_mm = permute(pixel_cp_mm{s}, [4,1,2,3]);
Y_pct = permute(pixel_cp_pct{s}, [4,1,2,3]);

% 数据扁平化
Y_mm = reshape(Y_mm, length(years), []);
Y_pct = reshape(Y_pct, length(years), []);

% 线性回归
beta_mm = X \ Y_mm;
beta_pct = X \ Y_pct;

% 恢复空间结构
pixel_kcp_mm{s} = reshape(beta_mm(2,:), lats, lons, 4);
pixel_kcp_pct{s} = reshape(beta_pct(2,:), lats, lons, 4);

% 季节处理完成清理
clear Y_* beta_*

end

%% 存储数据
save ('/workdir/tanjy/run/result_era5/pixel_cp.mat', 'pixel_cp_mm', 'pixel_cp_pct', '-v7.3')
save ('/workdir/tanjy/run/result_era5/pixel_kcp.mat', 'pixel_kcp_mm', 'pixel_kcp_pct', '-v7.3')

Thought

Clarifying precipitation calculations

The code processes ERA5 convective precipitation NetCDF files, calculating the precipitation amount for each grid point, broken down by local time bin for each season from 1980 to 2024. It then computes the fraction of precipitation in each time bin, calculates trends (slopes) for both the amount and percentage of precipitation over time, and saves those results for analysis. This provides insight into how convective precipitation has changed across seasons and years.

Thought for 6 seconds

Answer

这段 MATLAB 脚本的总体功能是：对 ERA5 再分析数据中的对流性降水（convective precipitation, “cp”）进行分季节、分地方时段的累积统计，计算其在各时段的占比，并对年际变化趋势（斜率）进行线性回归分析，最后将结果保存到 MAT 文件中。

具体流程可分为以下几部分：

  1. 基础参数与路径设置

    • 定义分析年份 years = 1980:2024，地方时段划分为 0–6、6–12、12–18、18–24 小时（time_bins = 0:6:24），以及 4 个季节（春、夏、秋、冬）和全年。
    • 指定数据根目录 data_root 和文件名格式 cp_%d%02d.nc。

  2. 地理网格与掩膜

    • 从本地 MAT 文件中加载陆地掩膜 landmask_era5 和海岸线数据。
    • 构建纬度（90 至 –60）和经度（0 至 359.75）向量，以及它们的长度 lats, lons。

  3. 初始化存储结构

    • 用 cell 数组 pixel_cp_mm{1:5} 分别存储春、夏、秋、冬、全年 4 个时段的降水累积量（单位 mm）。数组维度为 [纬度 &times; 经度 &times; 时段(4) &times; 年份数]。

  4. 主循环：读取 NetCDF 并累积降水

    • 遍历 years_data = 1980:2025，跳过 1980 年 1–2 月和 2025 年 3–12 月，以完整覆盖跨年冬季。
    • 对每个月：
      1. 用 ncread 读取时间（秒自 1970-01-01）并转换为 datetime，读取 cp 变量（原始维度 [lon,lat,time]，用 permute 调整为 [lat,lon,time]）。
      2. 对每个时刻 t：
        • 提取 UTC 小时 utc_hour = hour(dt)，并通过 local_time = mod(utc_hour + lon/15,24) 计算地方时。
        • 对 4 个时段分别做掩膜：只保留陆地像素且落在该时段内的降水量 precip_h。
        • 判断该时刻所属季节（MAM、JJA、SON、DJF），并根据月、年关系确定“季节年”（如 1、2 月属上一年冬季）。
        • 将 precip_h 加到对应季节、对应年份、对应时段的累积量上。

  5. 全年累积与百分比转换

    • 将前 4 个季节的降水累加到 cell{5}，得到全年降水累积。
    • 对每个季节，沿 4 个时段方向求和得到总降水 total，再用 bsxfun(@rdivide, …) 计算各时段占总降水的百分比 pixel_cp_pct。

  6. 年际趋势（斜率）分析

    • 构造回归设计矩阵 
      matlab
      X = [ones(length(years),1), years(:)-mean(years)];
      用中心化处理后的年份做线性回归。
    • 对每个季节：
      1. 将 pixel_cp_mm{s} 和 pixel_cp_pct{s} 从 [lat,lon,bin,year] 重排为 [year &times; (lat*lon*bin)] 的矩阵 Y_mm, Y_pct。
      2. 通过最小二乘 beta = X \ Y，提取第二行（斜率），再 reshape 回 [lat,lon,bin]，得到每个网格点、每个时段的年际变化率 pixel_kcp_mm（mm/年）和 pixel_kcp_pct（百分比点/年）。

  7. 结果保存

    • 将累积量与百分比结果保存到 
      matlab
      pixel_cp.mat    % 包含 pixel_cp_mm, pixel_cp_pct
pixel_kcp.mat   % 包含 pixel_kcp_mm, pixel_kcp_pct

简而言之，这段脚本是为了研究 1980–2024 年间 ERA5 对流性降水在不同季节和不同地方时段的分布特征，以及这种分布的长期趋势。

Share this Q&A