clear;
ncsu = [bign(5) bigc(7) bigs(9) bigu(11)];
newncsu = 20*ncsu;
[nx ny] = size(newncsu);
nx
ny
newncsu = newncsu(nx:-1:1,:);
newncsu1 = uint8(newncsu);
imwrite(newncsu1, 'ncsu.jpg');
u = newncsu;
for i = 1:nx
  for j = 1:ny
      % This is the big wave with periodic noise.
      %u(i,j) = 0+100*(1+sin(2*pi*((j-1)/ny)*5)) + ...
                    15.*(1+sin(2*pi*((i-1)/nx)*80))+...
                    15.*(1+sin(2*pi*((j-1)/ny)*80));
      % This is NCSU with periodic noise.
      u(i,j) = u(i,j) + ...
                    15.*(1+sin(2*pi*((i-1)/nx)*80))+...
                    15.*(1+sin(2*pi*((j-1)/ny)*80));
  end
end
sinencsu = uint8(u);
imwrite(sinencsu, 'sinencsu.jpg');
fftu = fft2(u,2*nx-1,2*ny-1);
fftu = fftshift(fftu);
subplot(2,2,1);
mesh(u');
subplot(2,2,2);
mesh(log(1+(abs(fftu))));
filter = ones(2*nx-1,2*ny-1);
d0 = 160;      % Use Butterworth band-reject filter.
n = 4;
w = 20;
for i = 1:2*nx-1
    for j =1:2*ny-1
        dist = ((i-(nx+1))^2 + (j-(ny+1))^2)^.5;
        if dist ~= d0
            filter(i,j)= 1/(1 + (dist*w/(dist^2 - d0^2))^(2*n));
        else
            filter(i,j) = 0;
        end
    end
end 
subplot(2,2,3);
mesh(filter);
fil_ncsu  = filter.*fftu;
subplot(2,2,4);
mesh(log(1+abs(fil_ncsu)));
fil_ncsu = ifftshift(fil_ncsu);
fil_ncsu = ifft2(fil_ncsu,2*nx-1,2*ny-1);
fil_ncsu = real(fil_ncsu(1:nx,1:ny));
fil_ncsu = uint8(fil_ncsu);
imwrite(fil_ncsu, 'sinencsu_fil.jpg');