Finding a points in 2 different matrices that satisfies a minimum criteria

Hi all,

Please refer to the .m files attached to this question - I would be grateful if you could actually try and run this code (and hopefully to come up with a solution) with the real data.

Need your assistance with this matter as I cant seem to figure out a rigorous solution to cover all cases (already asked this question here in the past with no avail).

Given the matrices r_sp, theta_sp and Dose_alpha_total. These matrices are pre-calcualted and given as is.

Main objective: finding the Dose_alpha_total value at a specific r & theta location, indicated by the user.

Catch 1: the values of r & theta that the user will enter, don't exist in r_sp & theta_sp matrices, so I need to find the cloesest value existing in both, and then performing a 2D interploation in order to give the appropriate dose value.

Catch 2: the nearest values found on r_sp & theta_sp must be on the same index on both matrices, as ultimately, the dose wil lbe calcualted by refering the this very same location on Dose_alpha_total matrix.

Catch 3: it is impossible to define a "tolerance" value to aid with the minimzation issue on Catch 1 above, as the differences between r_sp and theta_sp are dynamic, meaning that if one defines a tolerance X for one input-case which will work just fine, the very same X tolerance won't work with the next input-case, where the tolerance will shift by an order of magnitude.

This is the code I've been trying to implement, without any success so far:

Please note that the if-elseif-elseif-else bulk of code treats with finding a value larger/smaller than the one found on I_row1 & I_col1, so the interpolation would be possible to execute. I tried using scatteredInterpolant with no success.

The assignment of value 999 is to eliminate the cells that do not meet the if criteria, so the code wil be able to look for the next possible cell.

Also please note that while this code may work for certain values, it won't work for the case of, for instance, r_in=0.052 and theta_in = 0, as the r_diff & theta_diff matrices will find a cell that is extremely far from the input values, but its difference is indeed the lowest since the second parameter is extremely close, so it pulls the entire sum down.

I'm also quite sure that this code is extremely inefficient but it is what it is for nowr_diff=abs((r_sp- r_in)./r_sp);

theta_diff=abs((theta_sp -theta_in)./theta_sp);
add=r_diff+theta_diff;
[M,I] = min(add(:));
[I_row1, I_col1] = ind2sub(size(add),I);
dose1=Dose_alpha_total(I_row1,I_col1);
if r_sp(I_row1, I_col1)<r_in && theta_sp(I_row1, I_col1)< theta_in
    theta_sp(I_row1,I_col1)=999;
    r_sp(I_row1,I_col1)=999;
    r_diff=abs((r_sp- r_in)./r_sp);
    theta_diff=abs((theta_sp -theta_in)./theta_sp);
    add=r_diff+theta_diff;
    [M,I] = min(add(:));
    [I_row2, I_col2] = ind2sub(size(add),I);
    while r_sp(I_row2, I_col2)<r_in || theta_sp(I_row2, I_col2)< theta_in
        theta_sp(I_row2,I_col2)=999;
        r_sp(I_row2,I_col2)=999;
        r_diff=abs((r_sp- r_in)./r_sp);
        theta_diff=abs((theta_sp -theta_in)./theta_sp);
        add=r_diff+theta_diff;
        [M,I] = min(add(:));
        [I_row2, I_col2] = ind2sub(size(add),I);
    end
    dose2=Dose_alpha_total(I_row2,I_col2);
    r_sp=r_sp_f; %r_sp_f equals to the "untouched" r_sp matrix, so upon finding the wanted value, r_sp should return to its original, untouched conditions
    theta_sp=theta_sp_f; % same logic as r_sp_f
elseif r_sp(I_row1, I_col1)>r_in && theta_sp(I_row1, I_col1)> theta_in
        theta_sp(I_row1,I_col1)=999;
        r_sp(I_row1,I_col1)=999;
        r_diff=abs((r_sp- r_in)./r_sp);
        theta_diff=abs((theta_sp -theta_in)./theta_sp);
        add=r_diff+theta_diff;
        [M,I] = min(add(:));
        [I_row2, I_col2] = ind2sub(size(add),I);
        while r_sp(I_row2, I_col2)>r_in || theta_sp(I_row2, I_col2)> theta_in
            theta_sp(I_row2,I_col2)=999;
            r_sp(I_row2,I_col2)=999;
            r_diff=abs((r_sp- r_in)./r_sp);
            theta_diff=abs((theta_sp -theta_in)./theta_sp);
            add=r_diff+theta_diff;
            [M,I] = min(add(:));
            [I_row2, I_col2] = ind2sub(size(add),I);
        end 
        dose2=Dose_alpha_total(I_row2,I_col2);
        r_sp=r_sp_f;
        theta_sp=theta_sp_f;

elseif r_sp(I_row1, I_col1)<r_in && theta_sp(I_row1, I_col1)>theta_in
        theta_sp(I_row1,I_col1)=999;
        r_sp(I_row1,I_col1)=999;
        r_diff=abs((r_sp- r_in)./r_sp);
        theta_diff=abs((theta_sp -theta_in)./theta_sp);
        add=r_diff+theta_diff;
        [M,I] = min(add(:));
        [I_row2, I_col2] = ind2sub(size(add),I);
        while r_sp(I_row2, I_col2)<r_in || theta_sp(I_row2, I_col2)>theta_in
            theta_sp(I_row2,I_col2)=999;
            r_sp(I_row2,I_col2)=999;
            r_diff=abs((r_sp- r_in)./r_sp);
            theta_diff=abs((theta_sp -theta_in)./theta_sp);
            add=r_diff+theta_diff;
            [M,I] = min(add(:));
            [I_row2, I_col2] = ind2sub(size(add),I);
        end 
        dose2=Dose_alpha_total(I_row2,I_col2);
        r_sp=r_sp_f;
        theta_sp=theta_sp_f;
else 
    %r_sp(I_row1, I_col1)>r_in && theta_sp(I_row1, I_col1)< theta_in
    theta_sp(I_row1,I_col1)=999;
    r_sp(I_row1,I_col1)=999;
    r_diff=abs((r_sp- r_in)./r_sp);
    theta_diff=abs((theta_sp -theta_in)./theta_sp);
    add=r_diff+theta_diff;
    [M,I] = min(add(:));
    [I_row2, I_col2] = ind2sub(size(add),I);
    while r_sp(I_row2, I_col2)>r_in || theta_sp(I_row2, I_col2)< theta_in
        theta_sp(I_row2,I_col2)=999;
        r_sp(I_row2,I_col2)=999;
        r_diff=abs((r_sp- r_in)./r_sp);
        theta_diff=abs((theta_sp -theta_in)./theta_sp);
        add=r_diff+theta_diff;
        [M,I] = min(add(:));
        [I_row2, I_col2] = ind2sub(size(add),I);
    end
    dose2=Dose_alpha_total(I_row2,I_col2);
    r_sp=r_sp_f;
    theta_sp=theta_sp_f;
end

dose_avg=(dose1+dose2)/2; % just to see how far the avg is from the interpolated value

dose_interp1=(((r_sp_f(I_row2,I_col2)-r_in)*(theta_sp_f(I_row2,I_col2)-theta_in))/((r_sp_f(I_row2,I_col2)-r_sp_f(I_row1,I_col1))*(theta_sp_f(I_row2,I_col2)-theta_sp_f(I_row1,I_col1))))*Dose_alpha_total(I_row1,I_col1);
dose_interp2=(((r_in-r_sp_f(I_row1,I_col1))*(theta_sp_f(I_row2,I_col2)-theta_in))/((r_sp_f(I_row2,I_col2)-r_sp_f(I_row1,I_col1))*(theta_sp_f(I_row2,I_col2)-theta_sp_f(I_row1,I_col1))))*Dose_alpha_total(I_row2,I_col1);
dose_interp3=(((r_sp_f(I_row2,I_col2)-r_in)*(theta_in-theta_sp_f(I_row1,I_col1)))/((r_sp_f(I_row2,I_col2)-r_sp_f(I_row1,I_col1))*(theta_sp_f(I_row2,I_col2)-theta_sp_f(I_row1,I_col1))))*Dose_alpha_total(I_row1,I_col2);
dose_interp4=(((r_in-r_sp_f(I_row1,I_col1))*(theta_in-theta_sp_f(I_row1,I_col1)))/((r_sp_f(I_row2,I_col2)-r_sp_f(I_row1,I_col1))*(theta_sp_f(I_row2,I_col2)-theta_sp_f(I_row1,I_col1))))*Dose_alpha_total(I_row2,I_col2);
dose_interp=dose_interp1+dose_interp2+dose_interp3+dose_interp4;