%load Finv.txt
%load c.txt

Finv=makelut;  % genetrate LUT from known a,c, and s values
c=convert_Finv_to_c(Finv); % differentiate and interp to get certainty function

k=logspace(-3,0,31);  % linear exposure values from -30dB to 0dB
K=log(k);             % log exposure values

%% photoq cement
disp('determining size by loading first frame f1');
f1=loadframe(1);
[M,N]=size(f1);
qsum=zeros(M,N);
Csum=zeros(M,N);
for i=1:31
  disp(sprintf('loading frame %d',i));
  f_i = loadframe(i);
  %% in Matlab may need to convert uint8 to double (not necessary in Octave)
  q=exp(Finv(f_i+1))/k(i); % foreach i, get estimate q=q\hat_i with inverse F
  % in Matlab q will be same size as f_i, but in Octave need to reshape
  q=reshape(q,M,N);
  C=c(f_i+1); % certainty image for running total of certainties
  % in Matlab C will be same size as f_i, but in Octave need to reshape
  C=reshape(C,M,N);
  Csum = Csum + C; % Cumulative sum of Certainty images
  qsum=qsum+q.*C;
end
q=qsum./Csum; % final quantity of light is weighted sum

save cementedimage_with_known_a_c_and_s_parameters q