function HW1_Matlab_2020
close all  % Closes any open figure windows


%% Problem 7
disp('Solution to Problem 7: 501 equally spaced points from 0 to 2 pi')
t = linspace(0,2*pi,501)

% For those who prefer loops
% for i=1:501
%     t(i) = 2*pi*(i-1)/500;
% end

% [1 POINT]

%% Problem 8

% Create x and y
for i=1:length(t)
    if (t(i)<pi)
      x(i) = log(t(i)+2)*sin(t(i));
      y(i) = log(t(i)+2)*cos(t(i));
    else
      x(i) = log(t(i)+2-pi)*sin(t(i));
      y(i) = log(t(i)+2-pi)*cos(t(i)-pi);
    end

end

% Plot the graph
plot(x,y,'LineWidth',2,'Color',[1 0 0]);
axes1 = gca;
ylabel('y');
xlabel('x');
title('Plot for problem 8');
box(axes1,'on');
set(axes1,'FontSize',16,'XGrid','on','YGrid','on');

% [2 POINTS + 1 point credit for labeling the plot] 

%% Problem 9
figure
M = zeros(101);
M(25:75,25:75) = 1;
M(27:73,27:73) = 0.75;
imshow(M,'InitialMagnification',400)
title('Image for problem 9');

% [2 POINTS]

%% Problem 10

%
% Uncomment the lines below to display the MATLAB steganograph
% defImage = pow2(get(0,'DefaultImageCData'),47); 
% for shift = 10:63 
%     img = bitshift(defImage,shift);
%     imshow(img,'DisplayRange',[min(img(:)) max(img(:))],'InitialMagnification',400)
%     title(num2str(shift));
%     pause(1)
% end
% for i=1:101
%     for j=1:101
%         image(i,j) = 100*(j-1) + i-1;
%     end
% end
% figure
% imshow(image,'DisplayRange',[min(image(:)) max(image(:))],'InitialMagnification',400)
%

cover = imread('cover_image.jpg');
payload = imread('payload_image.jpg');

package = hide_image(cover,payload);

figure
imshow(package)
title('Package image')
payload = recover_image(package,cover);
figure
imshow(payload)
title('Recovered Payload image')

%[5 POINTS] - for correctly implementing functions (see below)


%% Problem 11

k1 = 1;
k2 = 1;
k3 = 1;
k4 = 1;
nA = 1;
nB = 3;

tstart = 0;
tstop = 50;

initial_w = [1;1];
[times,sols] = ode45(@(t,w) brusselator(t,w,nA,nB,k1,k2,k3,k4),[tstart,tstop],initial_w);

figure
  plot1 = plot(times,sols(:,1:2))
  set(plot1(1),'DisplayName','n_C');
  set(plot1(2),'DisplayName','n_D');
  axes1 = gca;
  ylabel('Concentrations (mols/l)');
  xlabel('Time (hrs)');
  title('Brusselator');
  box(axes1,'on');
  set(axes1,'FontSize',16,'XGrid','on','YGrid','on');
  legend(axes1,'show');

  % [4 POINTS]
  

end
function package = hide_image(cover,payload)
% Function to hide a message or image inside another image
% cover: the raw cover image
% payload: the message or image to be hidden
% package: the modified image that includes the payload
%
% The payload is hidden as follows;
%     (1) Convert the cover image to a double precision matrix
%     (2) Store the cover in the first 2 sig fig of package
%     (3) Store the payload in the remaining sig fig.
%
% Check % rows & columns of payload and cover
[nc,mc] = size(cover);
[np,mp] = size(payload);
% If the payload image is larger than the cover, we will scale the cover
% to be at least as big as the payload.  This is not needed for the sample
% images provided.
scale_factor = max([np/nc,mp/mc,1]);
cover = imresize(cover,scale_factor);
% Convert the images to double precision
cover = double(cover);
payload = double(payload);
% Scale cover image to range 0<cover<0.99
% You could do this with a loop too.   
package = 0.99*(cover - min(cover(:)))/(max(cover(:))-min(cover(:)));
% Scale payload image to range 0<payload<0.01
% You could also scale the image using
% payload = 0.01*cover/max(max(cover)); as suggested in the HW.
scaled_payload = 0.01*(payload - min(payload(:)))/(max(payload(:))-min(payload(:)));
% Add the payload to the top left corner of the package
package(1:np,1:mp) = package(1:np,1:mp) + scaled_payload;

% Alternative simpler code that only works on images with the same size:
% cover = double(cover); % Convert to DP as suggested in HW
% payload = double(payload); % Convert to DP as suggested in HW
% scaled_cover = cover/max(max(cover)); % as suggested in the HW.
% scaled_payload = payload/max(max(payload)); % Following suggestion in HW
% package = 0.99*scaled_cover + 0.01*scaled_payload;

end
function payload = recover_image(package,cover)
% Function to recover a message or image embedded inside another image
% cover: the raw cover image
% package: the modified image that includes the payload
% payload: the recovered message or image
% Check # rows and columns of the cover image and package; resize cover
[nc,mc] = size(cover);
[np,mp] = size(package);
scale_factor = max([min([np/nc,mp/mc]),1]);
cover = imresize(cover,scale_factor);
% Convert the images to double precision
cover = double(cover);
scaled_cover = 0.99*(cover - min(cover(:)))/(max(cover(:))-min(cover(:)));
% The payload is stored in the last sig figs of the top left corner of the cover
% To recover the payload we just subtract the cover image and
% scale appropriately.
payload = 100*(package-scaled_cover(1:np,1:mp));

% Alternative simpler code that only works on images with the same size:
% Convert the images to double precision
% cover = double(cover);
% scaled_cover = cover/max(max(cover)); % as suggested in the HW.
% payload = 100*(package - 0.99*scaled_cover);

end
function dwdt = brusselator(t,w,nA,nB,k1,k2,k3,k4)
% Function to calculate derivatives to solve Brusselator ODE

nC = w(1); nD = w(2);

dnCdt = k1*nA + k2*nC^2*nD - k3*nB*nC - k4*nC;
dnDdt = k3*nB*nC - k2*nC^2*nD;

dwdt = [dnCdt;dnDdt];

end

