nicklasorte · nicklasorte · Mar 27, 2026
diff --git a/benchmark_subchunk_agg_check_maxazi_rev10_rev13_real.m b/benchmark_subchunk_agg_check_maxazi_rev10_rev13_real.m
@@ -0,0 +1,46 @@
+function results = benchmark_subchunk_agg_check_maxazi_rev10_rev13_real(app,cell_aas_dist_data,array_bs_azi_data,radar_beamwidth,min_azimuth,max_azimuth,base_protection_pts,point_idx,on_list_bs,cell_sim_chunk_idx,rand_seed1,agg_check_reliability,on_full_Pr_dBm,clutter_loss,custom_antenna_pattern,sub_point_idx)
+%BENCHMARK_SUBCHUNK_AGG_CHECK_MAXAZI_REV10_REV13_REAL
+% Runtime benchmark for rev10 vs rev13 using identical real inputs.
+
+if exist('subchunk_agg_check_maxazi_rev10','file')~=2
+    error('benchmark_subchunk_agg_check_maxazi_rev10_rev13_real:MissingRev10', ...
+        'subchunk_agg_check_maxazi_rev10.m was not found on MATLAB path.');
+end
+if exist('subchunk_agg_check_maxazi_rev13','file')~=2
+    error('benchmark_subchunk_agg_check_maxazi_rev10_rev13_real:MissingRev13', ...
+        'subchunk_agg_check_maxazi_rev13.m was not found on MATLAB path.');
+end
+
+opts = struct();
+opts.AziChunkRev13 = 32;
+
+fprintf('\n=== REV10 vs REV13 REAL-INPUT BENCHMARK ===\n');
+fprintf('AZI_CHUNK rev13: %d\n',opts.AziChunkRev13);
+
+rev10_tic=tic;
+out_rev10=subchunk_agg_check_maxazi_rev10(app,cell_aas_dist_data,array_bs_azi_data, ...
+    radar_beamwidth,min_azimuth,max_azimuth,base_protection_pts,point_idx,on_list_bs, ...
+    cell_sim_chunk_idx,rand_seed1,agg_check_reliability,on_full_Pr_dBm,clutter_loss, ...
+    custom_antenna_pattern,sub_point_idx); %#ok<NASGU>
+runtime_rev10=toc(rev10_tic);
+
+rev13_tic=tic;
+out_rev13=subchunk_agg_check_maxazi_rev13(app,cell_aas_dist_data,array_bs_azi_data, ...
+    radar_beamwidth,min_azimuth,max_azimuth,base_protection_pts,point_idx,on_list_bs, ...
+    cell_sim_chunk_idx,rand_seed1,agg_check_reliability,on_full_Pr_dBm,clutter_loss, ...
+    custom_antenna_pattern,sub_point_idx,opts.AziChunkRev13); %#ok<NASGU>
+runtime_rev13=toc(rev13_tic);
+
+speedup = runtime_rev10 ./ runtime_rev13;
+
+fprintf('Runtime rev10: %.6f s\n',runtime_rev10);
+fprintf('Runtime rev13: %.6f s\n',runtime_rev13);
+fprintf('Speedup rev10/rev13: %.3fx\n',speedup);
+
+results = struct();
+results.runtime_rev10 = runtime_rev10;
+results.runtime_rev13 = runtime_rev13;
+results.speedup = speedup;
+results.azi_chunk_rev13 = opts.AziChunkRev13;
+
+end
diff --git a/subchunk_agg_check_maxazi_rev13.m b/subchunk_agg_check_maxazi_rev13.m
@@ -0,0 +1,121 @@
+function [sub_array_agg_check_mc_dBm]=subchunk_agg_check_maxazi_rev13(app,cell_aas_dist_data,array_bs_azi_data,radar_beamwidth,min_azimuth,max_azimuth,base_protection_pts,point_idx,on_list_bs,cell_sim_chunk_idx,rand_seed1,agg_check_reliability,on_full_Pr_dBm,clutter_loss,custom_antenna_pattern,sub_point_idx,varargin)
+%SUBCHUNK_AGG_CHECK_MAXAZI_REV13 Monte Carlo aggregate check with tunable azimuth chunking.
+% rev13 goals:
+%   1) remove per-iteration RNG reseeding overhead;
+%   2) keep azimuth chunking as a memory/performance tuning knob;
+%   3) preserve rev9/rev10 output contract (max aggregate dBm over sim azimuth).
+
+% Tuning knob: larger chunks can improve compute throughput but may increase peak memory.
+AZI_CHUNK_DEFAULT=32; % Selected from real-input chunk sweep benchmark (best runtime vs tested sizes).
+DEBUG_CHECKS=false;
+azi_chunk=AZI_CHUNK_DEFAULT;
+if ~isempty(varargin)
+    azi_chunk=varargin{1};
+end
+azi_chunk=max(1,round(azi_chunk));
+
+array_aas_dist_data=cell_aas_dist_data{2};
+aas_dist_azimuth=cell_aas_dist_data{1};
+mod_azi_diff_bs=array_bs_azi_data(:,4);
+
+% Off-axis EIRP lookup at BS-relative azimuth.
+nn_azi_idx=nearestpoint_app(app,mod_azi_diff_bs,aas_dist_azimuth);
+super_array_bs_eirp_dist=array_aas_dist_data(nn_azi_idx,:);
+
+% Simulation azimuth grid.
+[array_sim_azimuth,num_sim_azi]=calc_sim_azimuths_rev3_360_azimuths_app(app,radar_beamwidth,min_azimuth,max_azimuth);
+
+% BS->point azimuths.
+sim_pt=base_protection_pts(point_idx,:);
+bs_azimuth=azimuth(sim_pt(1),sim_pt(2),on_list_bs(:,1),on_list_bs(:,2));
+
+% MC iteration indices for this sub-point.
+sub_mc_idx=cell_sim_chunk_idx{sub_point_idx}; %#ok<NASGU>
+num_mc_idx=length(sub_mc_idx);
+num_bs=length(bs_azimuth);
+sub_array_agg_check_mc_dBm=NaN(num_mc_idx,1);
+
+% -------------------------------------------------------------------------
+% STEP 1: MC random pre-generation using a single RNG seeding call.
+% Draw in [rel_min, rel_max] for PR, EIRP, clutter random reliabilities.
+% -------------------------------------------------------------------------
+rel_min=min(agg_check_reliability);
+rel_max=max(agg_check_reliability);
+
+if rel_min==rel_max
+    rand_pr_all=repmat(rel_min,num_bs,num_mc_idx);
+    rand_eirp_all=rand_pr_all;
+    rand_clutter_all=rand_pr_all;
+else
+    rng(rand_seed1);
+    rel_span=(rel_max-rel_min);
+    rand_pr_all=rel_min+rel_span.*rand(num_bs,num_mc_idx);
+    rand_eirp_all=rel_min+rel_span.*rand(num_bs,num_mc_idx);
+    rand_clutter_all=rel_min+rel_span.*rand(num_bs,num_mc_idx);
+end
+
+% -------------------------------------------------------------------------
+% STEP 2: Precompute off-axis gain matrix once for all (bs,sim_azimuth).
+% Keep nearestpoint semantics stable.
+% -------------------------------------------------------------------------
+pat_az=mod(custom_antenna_pattern(:,1),360);
+pat_gain=custom_antenna_pattern(:,2);
+[pat_az_unique,ia_unique]=unique(pat_az,'stable');
+pat_gain_unique=pat_gain(ia_unique);
+
+off_axis_gain_matrix=NaN(num_bs,num_sim_azi);
+for azimuth_idx=1:1:num_sim_azi
+    sim_azimuth=array_sim_azimuth(azimuth_idx);
+    rel_az=mod(bs_azimuth-sim_azimuth,360);
+    ant_deg_idx=nearestpoint_app(app,rel_az,pat_az_unique);
+    off_axis_gain_matrix(:,azimuth_idx)=pat_gain_unique(ant_deg_idx);
+end
+
+% -------------------------------------------------------------------------
+% STEP 3: RNG-free MC pathloss terms for each MC realization.
+% -------------------------------------------------------------------------
+sort_monte_carlo_pr_dBm_all=NaN(num_bs,num_mc_idx);
+for loop_idx=1:1:num_mc_idx
+    pre_sort_monte_carlo_pr_dBm=monte_carlo_Pr_dBm_rev2_app(app,agg_check_reliability,on_full_Pr_dBm,rand_pr_all(:,loop_idx));
+    rand_norm_eirp=monte_carlo_super_bs_eirp_dist_rev5(app,super_array_bs_eirp_dist,agg_check_reliability,rand_eirp_all(:,loop_idx));
+    monte_carlo_clutter_loss=monte_carlo_clutter_rev3_app(app,agg_check_reliability,clutter_loss,rand_clutter_all(:,loop_idx));
+
+    sort_monte_carlo_pr_dBm_all(:,loop_idx)=pre_sort_monte_carlo_pr_dBm+rand_norm_eirp-monte_carlo_clutter_loss;
+end
+
+% -------------------------------------------------------------------------
+% STEP 4: Aggregate over BS in watts, convert back to dBm, then max over az.
+% -------------------------------------------------------------------------
+for loop_idx=1:1:num_mc_idx
+    base_mc=sort_monte_carlo_pr_dBm_all(:,loop_idx);
+    max_azi_agg=-Inf;
+
+    for azi_start=1:azi_chunk:num_sim_azi
+        azi_end=min(azi_start+azi_chunk-1,num_sim_azi);
+        chunk_gain=off_axis_gain_matrix(:,azi_start:azi_end);
+        sort_temp_mc_dBm=base_mc+chunk_gain;
+
+        if DEBUG_CHECKS
+            if any(isnan(sort_temp_mc_dBm),'all')
+                error('subchunk_agg_check_maxazi_rev13:NaNTempDbm','NaN detected in sort_temp_mc_dBm');
+            end
+        end
+
+        binary_sort_mc_watts=db2pow(sort_temp_mc_dBm)/1000;
+        if DEBUG_CHECKS
+            if any(isnan(binary_sort_mc_watts),'all')
+                error('subchunk_agg_check_maxazi_rev13:NaNWatt','NaN detected in binary_sort_mc_watts');
+            end
+        end
+
+        azimuth_agg_dBm_chunk=pow2db(sum(binary_sort_mc_watts,1,'omitnan')*1000);
+        chunk_max=max(azimuth_agg_dBm_chunk,[],'omitnan');
+        if chunk_max>max_azi_agg
+            max_azi_agg=chunk_max;
+        end
+    end
+
+    sub_array_agg_check_mc_dBm(loop_idx,1)=max_azi_agg;
+end
+
+end
diff --git a/validate_subchunk_agg_check_maxazi_rev11_rev13_statistical.m b/validate_subchunk_agg_check_maxazi_rev11_rev13_statistical.m
@@ -0,0 +1,161 @@
+function results = validate_subchunk_agg_check_maxazi_rev11_rev13_statistical(app,cell_aas_dist_data,array_bs_azi_data,radar_beamwidth,min_azimuth,max_azimuth,base_protection_pts,point_idx,on_list_bs,cell_sim_chunk_idx,rand_seed1,agg_check_reliability,on_full_Pr_dBm,clutter_loss,custom_antenna_pattern,sub_point_idx)
+%VALIDATE_SUBCHUNK_AGG_CHECK_MAXAZI_REV11_REV13_STATISTICAL
+% Statistical and runtime comparison for rev11 vs rev13 using identical real inputs.
+
+if exist('subchunk_agg_check_maxazi_rev11','file')~=2
+    error('validate_subchunk_agg_check_maxazi_rev11_rev13_statistical:MissingRev11', ...
+        'subchunk_agg_check_maxazi_rev11.m was not found on MATLAB path.');
+end
+if exist('subchunk_agg_check_maxazi_rev13','file')~=2
+    error('validate_subchunk_agg_check_maxazi_rev11_rev13_statistical:MissingRev13', ...
+        'subchunk_agg_check_maxazi_rev13.m was not found on MATLAB path.');
+end
+
+% Fixed thresholds; fail closed on exceedance.
+opts = struct();
+opts.AziChunkRev11 = 128;
+opts.AziChunkRev13 = 32;
+opts.AbsDiffThreshold_dB = 0.50;
+opts.RelDiffThreshold = 0.05;
+opts.EnableP999 = true;
+
+fprintf('\n=== REV11 vs REV13 STATISTICAL VALIDATION ===\n');
+fprintf('AZI_CHUNK rev11: %d\n',opts.AziChunkRev11);
+fprintf('AZI_CHUNK rev13: %d\n',opts.AziChunkRev13);
+
+rev11_tic=tic;
+out_rev11=subchunk_agg_check_maxazi_rev11(app,cell_aas_dist_data,array_bs_azi_data, ...
+    radar_beamwidth,min_azimuth,max_azimuth,base_protection_pts,point_idx,on_list_bs, ...
+    cell_sim_chunk_idx,rand_seed1,agg_check_reliability,on_full_Pr_dBm,clutter_loss, ...
+    custom_antenna_pattern,sub_point_idx,opts.AziChunkRev11);
+runtime_rev11=toc(rev11_tic);
+
+rev13_tic=tic;
+out_rev13=subchunk_agg_check_maxazi_rev13(app,cell_aas_dist_data,array_bs_azi_data, ...
+    radar_beamwidth,min_azimuth,max_azimuth,base_protection_pts,point_idx,on_list_bs, ...
+    cell_sim_chunk_idx,rand_seed1,agg_check_reliability,on_full_Pr_dBm,clutter_loss, ...
+    custom_antenna_pattern,sub_point_idx,opts.AziChunkRev13);
+runtime_rev13=toc(rev13_tic);
+
+speedup = runtime_rev11 ./ runtime_rev13;
+
+x11=out_rev11(:);
+x13=out_rev13(:);
+finite_mask=isfinite(x11) & isfinite(x13);
+x11=x11(finite_mask);
+x13=x13(finite_mask);
+
+if isempty(x11)
+    error('validate_subchunk_agg_check_maxazi_rev11_rev13_statistical:NoFiniteSamples', ...
+        'No finite paired samples available for statistical comparison.');
+end
+
+metrics={'mean','std','min','max','median','p90','p95','p99'};
+q=[0.90 0.95 0.99];
+
+s11.mean=mean(x11,'omitnan');
+s11.std=std(x11,0,'omitnan');
+s11.min=min(x11,[],'omitnan');
+s11.max=max(x11,[],'omitnan');
+s11.median=median(x11,'omitnan');
+q11=quantile(x11,q);
+s11.p90=q11(1); s11.p95=q11(2); s11.p99=q11(3);
+
+s13.mean=mean(x13,'omitnan');
+s13.std=std(x13,0,'omitnan');
+s13.min=min(x13,[],'omitnan');
+s13.max=max(x13,[],'omitnan');
+s13.median=median(x13,'omitnan');
+q13=quantile(x13,q);
+s13.p90=q13(1); s13.p95=q13(2); s13.p99=q13(3);
+
+if opts.EnableP999 && numel(x11)>=1000
+    metrics=[metrics {'p99_9'}]; %#ok<AGROW>
+    s11.p99_9=quantile(x11,0.999);
+    s13.p99_9=quantile(x13,0.999);
+end
+
+diff_table=struct();
+pass_flags=true(1,numel(metrics));
+for k=1:1:numel(metrics)
+    m=metrics{k};
+    v11=s11.(m);
+    v13=s13.(m);
+    abs_diff=abs(v13-v11);
+    rel_diff=abs_diff/max(abs(v11),eps);
+    allowed=max(opts.AbsDiffThreshold_dB,opts.RelDiffThreshold*max(abs(v11),1));
+
+    diff_table.(m).rev11=v11;
+    diff_table.(m).rev13=v13;
+    diff_table.(m).abs_diff=abs_diff;
+    diff_table.(m).rel_diff=rel_diff;
+    diff_table.(m).allowed_abs=allowed;
+    diff_table.(m).pass=abs_diff<=allowed;
+
+    pass_flags(k)=diff_table.(m).pass;
+end
+
+tail_fields=intersect({'p95','p99','p99_9'},metrics,'stable');
+tail_check=struct();
+tail_pass=true;
+for k=1:1:numel(tail_fields)
+    tf=tail_fields{k};
+    abs_diff=diff_table.(tf).abs_diff;
+    allowed=diff_table.(tf).allowed_abs;
+    tail_check.(tf).abs_diff=abs_diff;
+    tail_check.(tf).allowed_abs=allowed;
+    tail_check.(tf).pass=abs_diff<=allowed;
+    tail_pass=tail_pass && tail_check.(tf).pass;
+end
+
+overall_pass=all(pass_flags) && tail_pass;
+
+fprintf('Runtime rev11: %.6f s\n',runtime_rev11);
+fprintf('Runtime rev13: %.6f s\n',runtime_rev13);
+fprintf('Speedup rev11/rev13: %.3fx\n',speedup);
+
+fprintf('\nMetric comparison (rev13 - rev11):\n');
+for k=1:1:numel(metrics)
+    m=metrics{k};
+    fprintf('  %-7s | rev11=%10.4f | rev13=%10.4f | abs=%.4f | allow=%.4f | %s\n', ...
+        m,diff_table.(m).rev11,diff_table.(m).rev13,diff_table.(m).abs_diff, ...
+        diff_table.(m).allowed_abs,passfail(diff_table.(m).pass));
+end
+
+fprintf('\nUpper-tail checks:\n');
+for k=1:1:numel(tail_fields)
+    tf=tail_fields{k};
+    fprintf('  %-7s | abs=%.4f | allow=%.4f | %s\n',tf,tail_check.(tf).abs_diff, ...
+        tail_check.(tf).allowed_abs,passfail(tail_check.(tf).pass));
+end
+
+if overall_pass
+    fprintf('\nPASS: rev13 is statistically equivalent to rev11 under configured thresholds.\n');
+else
+    fprintf('\nFAIL: rev13 drift exceeded configured thresholds.\n');
+    error('validate_subchunk_agg_check_maxazi_rev11_rev13_statistical:DriftExceeded', ...
+        'Fail-closed: statistical drift exceeded configured thresholds.');
+end
+
+results=struct();
+results.runtime_rev11_s=runtime_rev11;
+results.runtime_rev13_s=runtime_rev13;
+results.speedup_rev11_over_rev13=speedup;
+results.n_samples=numel(x11);
+results.metrics=metrics;
+results.summary_rev11=s11;
+results.summary_rev13=s13;
+results.diffs=diff_table;
+results.upper_tail=tail_check;
+results.thresholds=opts;
+results.pass=overall_pass;
+
+end
+
+function txt=passfail(tf)
+if tf
+    txt='PASS';
+else
+    txt='FAIL';
+end
+end