% plotf(a_vector_in, c_vector_in)
%   given two column vectors of a and c values,
%   plots the response functions on the same graph
%   this is a smooth curve that fits through the response function,
%   regluarizing (smoothing it),( remember that a=sqrt(e^a) and c=sqrt(c)
%   are the actual values being optimized, and should be passed in to this 
%   function: i.e. in this function, it calcs c=c*c, a=log(a*a), so a,c are as
%   solved by the comparasolve)

function plotf(a_vector_in, c_vector_in)
  size_a = size(a_vector_in);
  num_ac = size_a(1); 
  num_ac
  disp('for now, setting K0 to have 0 at maximal Q value')
  for i=1:num_ac
    a=a_vector_in(i,1);
    a=log(a*a)
    c=c_vector_in(i,1);
    c=c*c
    b=0; % without loss of generality
    K=1; % in this case, k=2 so K=1

    K0=-1

    Q=linspace(-20,10,1000); % Quantity (log quantity) of light on image sensor
    f(i,:)=1../((1+exp(-(a*(Q+K0)+b))).^c);
  end 

  %%plot(Q,256*f)
  size(f)

plot(Q( 1:16:1000),256*f(1, 1:16:1000),"-@14;k=1.50;",
     Q( 2:16:1000),256*f(2, 1:16:1000),"-;k=1.75;",
     Q(10:16:1000),256*f(3,10:16:1000),"-;k=2.00;",
     Q( 5:16:1000),256*f(4, 5:16:1000),"-;k=2.25;",
     Q( 1:16:1000),256*f(5, 1:16:1000),"-;k=2.50;",
     Q(14:16:1000),256*f(6,14:16:1000),"-;k=2.75;",
     Q(11:16:1000),256*f(7,11:16:1000),"-;k=3.00;",
     Q( 1:16:1000),256*f(8, 1:16:1000),"-;k=3.25;"
    )

  xlabel('Quantity of light, Q = log(q)')
  ylabel('Pixel integer') % since F^-1(Q)=1/(1+e^-Q)
  grid
endfunction
