fid=fopen('lena_rumorosa.raw','rb');
if fid<0
    error('Errore apertura');
end
x=(fread(fid,[256 256],'float'))';
figure;
imshow(x,[]);

k=15;
cu=1/(sqrt(2)*2);

sigmal=colfilt(x,[k k],'sliding',@std);
mul=colfilt(x,[k k],'sliding',@mean);

cx=sigmal./mul;
r=cu^2./(cx.^2);
wl=(1-r)./(1+cu^2);

if r>1
    y=mul;
else
    y=(x.*wl)+mul.*(1-wl);
end

imshow(y,[]);