% <MAKEFFCORR> - Creates a flatfield correction file
%
%   [Im, <Sigma>] = makeffcorr (ImStack, <MaskImg>)
%   Im     : Flatfield correction image, normalized around 1 to correct
%            a pixel's different quantum efficiency for incoming photons. 
%   Sigma  : Standard deviation of each pixel of Im (ideally, this reflects
%            the counting statistics of each pixel, normalized to the mean
%            of Im).
%   ImStack: 3D stack of all measured flatfield images.
%   MaskImg: Binary mask to select the pixels on which makeffcorr.m acts.
%            Note: Use e.g. histmask.m to generate such a mask.
%   <argument> depicts an optional argument.

%==========================================================================
% 
% FUNCTION: makeffcorr.m 
%           ============
%
% $Date: 2006/12/11 12:39:13 $
% $Author: pauli_s $
% $Revision: 1.8 $
% $Source: /import/cvs/X/PILATUS/App/lib/X_PILATUS_Matlab/makeffcorr.m,v $
% $Tag: $
%
%
% <MAKEFFCORR> - Creates a flatfield correction file
% 
% Author(s):            R. Herger (RH)
% Co-author(s):         C.M. Schlepuetz (CS)
%                       P.R. Willmott (PRW) 
% Address:              Swiss Light Source (SLS)
%                       Paul Scherrer Institute
%                       CH - 5232 Villigen PSI
% Created:              2005/11/09
% 
% Change Log:
% -----------
% 
% 2005/11/09 (RH):
% - start
%
% 2005/11/11 (RH):
% - finished first version
%
% 2005/11/14 (RH):
% - changed mean to median for ImWork
%
% 2006/02/14 (RH & DM & PRW): 
% - corrected algorithm for finding the mean in-plane intensity 


%==========================================================================
%  Main function - makeffcorr
%                  ==========

function [Im, Sigma] = makeffcorr (ImStack, MaskImg)

%----------------------
% check input arguments

% checkout dimensions of stack, create a vector dim containing the
% in-plane dimensions for convenience

dim = size(ImStack);
xdim = dim(2);
ydim = dim(1);
clear dim

% are there 1 to 2 input arguments?
error(nargchk(1, 2, nargin));

if (nargin < 2)
  MaskImg = ones([ydim xdim]);  % set maskimg to ones if mask not given 
end;

% check if ImStack has 3 dimensions
if (ndims(ImStack) ~= 3)
  error(strcat('Invalid input for ''Stack'' in function ', ...
               ' immerge.\nUse ''help immerge'' for further', ...
               ' information.\n'), ...
               '')
end;

%----------------------
% check output argument


% are  1 to 2 output arguments specified?
error(nargoutchk(1, 2, nargout));

%-----------
% makeffcorr

% Initialize some variables
MaskZero = zeros ([ydim xdim]);

% create an image and a sigma matrix that contain the means and the
% standard deviations along the third dimension of the stack, respectively.
ImWork = mean (ImStack, 3) .* MaskImg;
SigmaWork = std (ImStack, 0, 3) .* MaskImg;


% do the normalization to the mean in-plane of ImWork and SigmaWork 
% (omitting the zeros), create flatfield correction and Sigma matrix
indices = (ImWork == 0);
MaskZero (indices) = 1;     % sets positions of dead pixels to unity in MaskZero 
MaskNonzero = ~MaskZero;    % inverse of MaskZero 
inPlaneMean = mean(nonzeros(ImWork));
InvIm = ImWork / inPlaneMean;   % divide mean array by its in-plane mean 
% Now take the reciprocal of InvIm as our correction image. Because there
% are dead pixels, we need to add unity to these pixel intensities to avoid
% dividing by zero. After inversion, we anyway cut out these pixels by
% multiplying by MaskNonZero. 
Im = (1 ./ (InvIm + MaskZero)) .* MaskNonzero;
Sigma = SigmaWork/inPlaneMean;
% clear some variables
clear ImWork SigmaWork inPlaneMean xidm ydim;
clear MaskZero MaskNonzero

%==========================================================================
%
%---------------------------------------------------%
% emacs setup:  force text mode to get no help with %
%               indentation and force use of spaces %
%               when tabbing.                       %
% Local Variables:                                  %
% mode:text                                         %
% indent-tabs-mode:nil                              %
% End:                                              %
%---------------------------------------------------%
%
% $Log: makeffcorr.m,v $
% Revision 1.8  2006/12/11 12:39:13  pauli_s
% added tag
%
% Revision 1.7  2006/11/30 13:57:57  pauli_s
% Remove unneeded files
%
% Revision 1.6  2006/02/17 16:36:06  willmott
% Improved comments
%
% Revision 1.5  2006/02/17 16:01:35  herger
% cleaned up some debug statements
%
% Revision 1.4  2006/02/14 16:45:32  willmott
% corrected algorithm for finding mean in-plane intensity
%
% Revision 1.3  2005/11/16 09:47:01  herger
% removed bug in mean value calculation (of variables ImWork and SigmaWork)
%
% Revision 1.2  2005/11/14 13:32:35  herger
% changed mean to median for ImWork
%
%
%================================= End of makeffcorr.m ====================