%OUT = COMPUTE1(OPERATION,IN)
%    This is the basis function for all the overloaded operators and
%    functions.
%OUT = COMPUTE1(OPERATION,IN,DIMS,OUT_DIP_TYPE,OUT_PHYS)
%    This is the form if you already called DO1INPUT(). Use this if
%    the format returned is not correct.

% (C) Copyright 1999-2009               Pattern Recognition Group
%     All rights reserved               Faculty of Applied Physics
%                                       Delft University of Technology
%                                       Lorentzweg 1
%                                       2628 CJ Delft
%                                       The Netherlands
%
% Cris Luengo, April 2000.
% 1 June 2000:     Using MATTYPE created this week.
% 2 June 2000:     Using FEVAL instead of BUILTIN.
% 18 January 2001: Added DIMS parameter. Changed all calls to this
%                  function accordingly.
% 14 August 2006:  Fixed output image being reshaped if input is 1xN.
% February 2008:   Adding pixel dimensions and units to dip_image. (BR)
% 11 March 2008:   Since MATLAB 7.0, computations can be done with singles.
% 17 December 2009: Using new function MATLABVER_GE.

function out = compute1(operation,in,dims,out_dip_type,out_phys)
if nargin == 2
   [in,dims,out_dip_type,out_phys] = do1input(in);
elseif nargin ~= 5
   error('Erroneus input.')
end
[out_type,complexity] = di_mattype(out_dip_type);
direct = 0;
if isa(in,'double')
   direct = 1;
else
   dosingle = 0;
   if matlabver_ge([7,0])
      if ~strcmp(out_type,'double')
         dosingle = 1;
      end
   end
   if dosingle & isa(in,'single')
      direct = 1;
   else
      s = size(in);
      N = prod(s);
      buffersize = dipgetpref('ComputationLimit');
      if dosingle
         buffersize = floor(buffersize/4);
         if N <= buffersize
            direct = 1;
            in = single(in);
         end
      else
         buffersize = floor(buffersize/8);
         if N <= buffersize
            direct = 1;
            in = double(in);
         end
      end
   end
end
if direct
   out = feval(operation,in);
   %#function int8, uint8, int16, uint16, int32, uint32, int64, uint64, single, double
   out = feval(out_type,out);
else
   out = di_create(s,out_type);
   if complexity, out = complex(out); end
   s = 1:buffersize:N;
   %disp(['We''re doing this in ',num2str(length(s)),' parts!']);
   if dosingle
      for ii=1:length(s)-1
         out(s(ii):s(ii+1)-1) = feval(operation,single(in(s(ii):s(ii+1)-1)));
      end
      out(s(end):N) = feval(operation,single(in(s(end):N)));
   else
      for ii=1:length(s)-1
         out(s(ii):s(ii+1)-1) = feval(operation,double(in(s(ii):s(ii+1)-1)));
      end
      out(s(end):N) = feval(operation,double(in(s(end):N)));
   end
end
% The next step not only creates the dip_image object, but also sets the
% datatype correctly. That is, we expect no conversion is done anymore,
% but a binary image will be labelled as such.
out = dip_image('trust_me',out,out_dip_type,dims,out_phys);
