Bug in phase screen generation when seed_density not equal dx_pixel

[nobias]

Hi! Thanks for sharing this nice approach and code! I think I found a bug in RandPhaseScreen_RealSpaceNM.m for the case when seed_density is not equal to dx_pixel.

Your current code, starting in line 52, is:

if (Nsame>1)
   for ii=1:Nsame:N(1)-Nsame
       for jj=1:Nsame:N(2)-Nsame
          ph_mask(ii+1:ii+Nsame,jj+1:jj+Nsame) =ph_mask(ii,jj);
       end
   end 
end

This generates some diagonal pattern, which probably is because you overwrite ph_mask with itself while you are changing it. I found that the following code fixes the issue:

ph_mask_same = zeros(size(ph_mask));
if Nsame > 1
    for ii = 1 : Nsame : N(1)-Nsame
        for jj = 1 : Nsame : N(2)-Nsame
            ph_mask_same(ii+1:ii+Nsame, jj+1:jj+Nsame) = ph_mask(ii,jj);            
        end
    end
    ph_mask = ph_mask_same;
end

Perhaps you would like to update the file.

[cxjphy]

Hi! Thanks for sharing this nice approach and code! I think I found a bug in RandPhaseScreen_RealSpaceNM.m for the case when seed_density is not equal to dx_pixel.

Your current code, starting in line 52, is:

if (Nsame>1)
   for ii=1:Nsame:N(1)-Nsame
       for jj=1:Nsame:N(2)-Nsame
          ph_mask(ii+1:ii+Nsame,jj+1:jj+Nsame) =ph_mask(ii,jj);
       end
   end 
end

This generates some diagonal pattern, which probably is because you overwrite ph_mask with itself while you are changing it. I found that the following code fixes the issue:

ph_mask_same = zeros(size(ph_mask));
if Nsame > 1
    for ii = 1 : Nsame : N(1)-Nsame
        for jj = 1 : Nsame : N(2)-Nsame
            ph_mask_same(ii+1:ii+Nsame, jj+1:jj+Nsame) = ph_mask(ii,jj);            
        end
    end
    ph_mask = ph_mask_same;
end

Perhaps you would like to update the file.

Thanks! We will take a look