view m-toolbox/classes/+utils/@math/intfact.m @ 33:5e7477b94d94 database-connection-manager

Add known repositories list to LTPDAPreferences
author Daniele Nicolodi <nicolodi@science.unitn.it>
date Mon, 05 Dec 2011 16:20:06 +0100
parents f0afece42f48
children
line wrap: on
line source

% INTFACT computes integer factorisation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% INTFACT tries to find two integers, P and Q, that satisfy
%
%          y = P/Q * x
%
%   >> [p,q] = intfact(y,x)
%
% The following call returns a parameter list object that contains the
% default parameter values:
%
% >> pl = intfact(utils, 'Params')
%
% The following call returns a string that contains the routine CVS version:
%
% >> version = intfact(utils,'Version')
%
% The following call returns a string that contains the routine category:
%
% >> category = intfact(utils,'Category')
%
% VERSION: $Id: intfact.m,v 1.9 2009/02/06 12:03:52 hewitson Exp $
%
% M Hewitson 26-05-08
%
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function varargout = intfact(varargin)
  
  % Get input sample rates
  fs2 = floor(1e10*varargin{1});
  fs1 = floor(1e10*varargin{2});
  
  g = gcd(fs2,fs1);
  
  varargout{1} = fs2/g;
  varargout{2} = fs1/g;
  
  return
  
  % Get min and max fs
  x = max(fs1,fs2);
  y = min(fs1,fs2);
  
  % If we already have two integers then it's easy
  if rem(fs2,1) == 0 && rem(fs1,1) == 0
    P = fs2;
    Q = fs1;
    % Get the greatest common divisor
    g = gcd(P,Q);
    % factor that out
    P = P / g;
    Q = Q / g;
    % If fs1 and fs2 are factors of each other, then it's easy
  elseif  rem(fs1,fs2) == 0
    P = fs1/fs2;
    Q = 1;
  elseif rem(fs2,fs1) == 0
    Q = fs2/fs1;
    P = 1;
    % Otherwise we search for two integers
  else
    % Start from N=0
    N = 0;
    s = x * 10^N;
    % Now compute initial P and Q
    if fs2 > fs1
      P  = 10^N;
      Q  = s / y;
    else
      Q  = 10^N;
      P  = s / y;
    end
    % Now look for the integers
    while rem(s,1) > 0 || rem(P,1) > 0 || rem(Q,1) > 0
      s = x * 10^N;
      % Now compute P and Q
      if fs2 > fs1
        P  = 10^N;
        Q  = s / y;
      else
        Q  = 10^N;
        P  = s / y;
      end
      N = N+1;
    end
    
    % Get the greatest common divisor
    g = gcd(P,Q);
    
    % factor that out
    P = P / g;
    Q = Q / g;
  end
  
  % Set outputs
  if nargout == 2
    varargout{1} = P;
    varargout{2} = Q;
  else
    error('### Incorrect number of outputs');
  end
  
end
% END