Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
56 changes: 56 additions & 0 deletions monte_carlo_super_bs_eirp_dist_rev7.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,56 @@
function [rand_norm_eirp]=monte_carlo_super_bs_eirp_dist_rev7(app,super_array_bs_eirp_dist,reliability,rand_numbers)
%MONTE_CARLO_SUPER_BS_EIRP_DIST_REV7 Correctness-first RNG-free MC EIRP interpolation.
% rev7 policy:
% - preserve rev5 interpolation semantics first;
% - keep shape-safe guards to prevent silent broadcasting/shape drift;
% - avoid aggressive spline PP rewrites unless helper validator proves exactness.

% Keep signature compatibility with prior revisions.
% app is intentionally unused.

DEBUG_CHECKS=false;

[num_rows,num_cols]=size(super_array_bs_eirp_dist);

if num_cols<=1
rand_norm_eirp=zeros(num_rows,1);
return;
end

rel_row=reliability(:).';
if numel(rel_row)~=num_cols
error('monte_carlo_super_bs_eirp_dist_rev7:ReliabilityLengthMismatch', ...
'reliability length (%d) must match number of EIRP columns (%d).',numel(rel_row),num_cols);
end

if ~issorted(rel_row)
[rel_row,sort_idx]=sort(rel_row,'ascend');
super_array_bs_eirp_dist=super_array_bs_eirp_dist(:,sort_idx);
end

xi=rand_numbers(:);
if numel(xi)~=num_rows
error('monte_carlo_super_bs_eirp_dist_rev7:QueryLengthMismatch', ...
'rand_numbers length (%d) must match number of EIRP rows (%d).',numel(xi),num_rows);
end

rel_min=rel_row(1);
rel_max=rel_row(end);
xi=min(max(xi,rel_min),rel_max);

if DEBUG_CHECKS
assert(isrow(rel_row),'Expected reliability to be a row vector for interp1.');
assert(size(super_array_bs_eirp_dist,1)==num_rows,'Y row size changed unexpectedly.');
assert(size(super_array_bs_eirp_dist,2)==numel(rel_row),'Y columns must match X length.');
if any(~isfinite(rel_row)) || any(~isfinite(super_array_bs_eirp_dist),'all') || any(~isfinite(xi))
error('monte_carlo_super_bs_eirp_dist_rev7:NonFiniteInput', ...
'Non-finite values detected in interpolation inputs.');
end
end

rand_norm_eirp=NaN(num_rows,1);
for n=1:1:num_rows
rand_norm_eirp(n)=interp1(rel_row,super_array_bs_eirp_dist(n,:),xi(n),'spline');
end

end
120 changes: 120 additions & 0 deletions subchunk_agg_check_maxazi_rev15.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,120 @@
function [sub_array_agg_check_mc_dBm]=subchunk_agg_check_maxazi_rev15(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_REV15 Monte Carlo aggregate check (rev11-compatible).
% rev15 policy:
% - keep rev11 RNG strategy/chunking/aggregation unchanged;
% - replace only EIRP helper call with correctness-first rev7 helper.

% Tuning knob: larger chunks can improve compute throughput but may increase peak memory.
AZI_CHUNK_DEFAULT=128;
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_rev7(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_rev15: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_rev15: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
184 changes: 184 additions & 0 deletions validate_monte_carlo_super_bs_eirp_dist_rev5_rev6.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,184 @@
function results=validate_monte_carlo_super_bs_eirp_dist_rev5_rev6(app,cell_aas_dist_data,array_bs_azi_data,agg_check_reliability,rand_seed1,cell_sim_chunk_idx,sub_point_idx,varargin)
%VALIDATE_MONTE_CARLO_SUPER_BS_EIRP_DIST_REV5_REV6
% Helper-level comparator for rev5 vs rev6 using identical helper-call inputs.
% Fail-closed policy: throw if helper drift exceeds tight thresholds.

must_exist('monte_carlo_super_bs_eirp_dist_rev5','MissingRev5');
must_exist('monte_carlo_super_bs_eirp_dist_rev6','MissingRev6');

opts=struct();
opts.NumSamples=[]; % default: infer from cell_sim_chunk_idx{sub_point_idx}
opts.MaxAbsDiffTol=1e-9;
opts.MeanAbsDiffTol=1e-10;
opts.RelDiffTol=1e-9;
opts.ReportTopK=5;
if ~isempty(varargin)
user_opts=varargin{1};
if isstruct(user_opts)
opts=merge_struct(opts,user_opts);
end
end

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);
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,:);

num_bs=size(super_array_bs_eirp_dist,1);
if isempty(opts.NumSamples)
sub_mc_idx=cell_sim_chunk_idx{sub_point_idx};
num_samples=length(sub_mc_idx);
else
num_samples=max(1,round(opts.NumSamples));
end

rel_min=min(agg_check_reliability);
rel_max=max(agg_check_reliability);
if rel_min==rel_max
rand_eirp_all=repmat(rel_min,num_bs,num_samples);
else
rng(rand_seed1);
rand_eirp_all=rel_min+(rel_max-rel_min).*rand(num_bs,num_samples);
end

x=agg_check_reliability(:).';
y=super_array_bs_eirp_dist;
query=rand_eirp_all;

fprintf('\n=== HELPER VALIDATION: rev5 vs rev6 ===\n');
fprintf('X (reliability grid) shape: [%d x %d]\n',size(x,1),size(x,2));
fprintf('Y (EIRP values) shape: [%d x %d]\n',size(y,1),size(y,2));
fprintf('Query shape: [%d x %d]\n',size(query,1),size(query,2));

out5=NaN(num_bs,num_samples);
out6=NaN(num_bs,num_samples);
for k=1:1:num_samples
out5(:,k)=monte_carlo_super_bs_eirp_dist_rev5(app,y,x,query(:,k));
out6(:,k)=monte_carlo_super_bs_eirp_dist_rev6(app,y,x,query(:,k));
end

delta=out6-out5;
abs_delta=abs(delta);

max_abs_diff=max(abs_delta,[],'all');
mean_abs_diff=mean(abs_delta,'all','omitnan');
den=max(abs(out5),1e-12);
rel_delta=abs_delta./den;
max_rel_diff=max(rel_delta,[],'all');
mean_rel_diff=mean(rel_delta,'all','omitnan');

lin_abs=abs_delta(:);
[sorted_abs,sorted_idx]=sort(lin_abs,'descend'); %#ok<ASGLU>
report_k=min(opts.ReportTopK,numel(sorted_idx));
worst=struct('linear_idx',cell(report_k,1),'row',[],'col',[], ...
'query_value',[],'rev5',[],'rev6',[],'abs_diff',[],'rel_diff',[]);
for i=1:report_k
idx=sorted_idx(i);
[r,c]=ind2sub(size(abs_delta),idx);
worst(i).linear_idx=idx;
worst(i).row=r;
worst(i).col=c;
worst(i).query_value=query(r,c);
worst(i).rev5=out5(r,c);
worst(i).rev6=out6(r,c);
worst(i).abs_diff=abs_delta(r,c);
worst(i).rel_diff=rel_delta(r,c);
end

x_sorted=sort(x,'ascend');
x_unique=unique(x_sorted,'stable');
x_min=x_sorted(1);
x_max=x_sorted(end);
if numel(x_unique)>=3
typical_step=median(diff(x_unique));
else
typical_step=max(x_max-x_min,eps);
end
endpoint_tol=max(typical_step*0.5,1e-12);
breakpoint_tol=max(typical_step*0.05,1e-12);
clamp_tol=1e-12;

near_endpoint=(abs(query-x_min)<=endpoint_tol) | (abs(query-x_max)<=endpoint_tol);
near_breakpoint=is_near_breakpoint(query,x_unique,breakpoint_tol);
is_clamped=(abs(query-rel_min)<=clamp_tol) | (abs(query-rel_max)<=clamp_tol);

sig_mask=abs_delta>max(opts.MaxAbsDiffTol,opts.RelDiffTol*max(abs(out5),1));
cluster=struct();
cluster.total_sig=nnz(sig_mask);
if cluster.total_sig>0
cluster.endpoint_fraction=nnz(sig_mask & near_endpoint)/cluster.total_sig;
cluster.breakpoint_fraction=nnz(sig_mask & near_breakpoint)/cluster.total_sig;
cluster.clamped_fraction=nnz(sig_mask & is_clamped)/cluster.total_sig;
else
cluster.endpoint_fraction=0;
cluster.breakpoint_fraction=0;
cluster.clamped_fraction=0;
end

fprintf('max abs diff: %.6e\n',max_abs_diff);
fprintf('mean abs diff: %.6e\n',mean_abs_diff);
fprintf('max rel diff: %.6e\n',max_rel_diff);
fprintf('mean rel diff: %.6e\n',mean_rel_diff);
fprintf('Significant-diff clustering: endpoints=%.1f%%, breakpoints=%.1f%%, clamped=%.1f%%\n', ...
100*cluster.endpoint_fraction,100*cluster.breakpoint_fraction,100*cluster.clamped_fraction);

for i=1:report_k
w=worst(i);
fprintf(' worst[%d] row=%d col=%d query=%.6g rev5=%.6g rev6=%.6g abs=%.3e rel=%.3e\n', ...
i,w.row,w.col,w.query_value,w.rev5,w.rev6,w.abs_diff,w.rel_diff);
end

pass=(max_abs_diff<=opts.MaxAbsDiffTol) && (mean_abs_diff<=opts.MeanAbsDiffTol) && (max_rel_diff<=opts.RelDiffTol);
fprintf('Helper comparison result: %s\n',passfail(pass));

results=struct();
results.options=opts;
results.input_shapes=struct('X',size(x),'Y',size(y),'Query',size(query));
results.thresholds=struct('max_abs',opts.MaxAbsDiffTol,'mean_abs',opts.MeanAbsDiffTol,'max_rel',opts.RelDiffTol);
results.max_abs_diff=max_abs_diff;
results.mean_abs_diff=mean_abs_diff;
results.max_rel_diff=max_rel_diff;
results.mean_rel_diff=mean_rel_diff;
results.worst_cases=worst;
results.cluster=cluster;
results.pass=pass;

if ~pass
error('validate_monte_carlo_super_bs_eirp_dist_rev5_rev6:DriftExceeded', ...
['Fail-closed: helper drift exceeded thresholds. ' ...
'max_abs=%.3e (tol %.3e), mean_abs=%.3e (tol %.3e), max_rel=%.3e (tol %.3e).'], ...
max_abs_diff,opts.MaxAbsDiffTol,mean_abs_diff,opts.MeanAbsDiffTol,max_rel_diff,opts.RelDiffTol);
end

end

function out=is_near_breakpoint(query,x_unique,tol)
out=false(size(query));
for i=1:1:numel(x_unique)
out=out | abs(query-x_unique(i))<=tol;
end
end

function must_exist(fname,errid)
if exist(fname,'file')~=2
error(['validate_monte_carlo_super_bs_eirp_dist_rev5_rev6:' errid], ...
'%s.m was not found on MATLAB path.',fname);
end
end

function out=merge_struct(base,override)
out=base;
fn=fieldnames(override);
for i=1:numel(fn)
out.(fn{i})=override.(fn{i});
end
end

function txt=passfail(tf)
if tf
txt='PASS';
else
txt='FAIL';
end
end
Loading