%MAKEWEIGHT  1st order system with given DC gain, crossover frequency,
%   and high frequency gain.
%
%   G = MAKEWEIGHT(DC,CROSSW,HF) creates a stable, 1st-order continuous
%   time state-space system G.  The frequency response of G satisfies
%   G(j*0) = DC, |G(j*CROSSW)| = 1, and G(j*infty) = HF.  It must be
%   that |DC| < 1 < |HF|, or |HF| < 1 < |DC|.
%
%   G = MAKEWEIGHT(DC,CROSSW,HF,TS) creates a stable, 1st-order discrete
%   time state-space system G, with sample-time TS.  The frequency
%   response of G satisfies G(exp(j*0*TS)) = DC, |G(exp(j*CROSSW*TS))| = 1,
%   and G(exp(j*PI)) = HF.  CROSSW*TS must be less than PI.
%
%   In both continuous-time and discrete-time, the low frequency
%   and high frequency targets must satisfy
%   |DC| < 1 < |HF|, or |HF| < 1 < |DC|.
%

function G = makeweight(dcval,bw,ginf,TS)

if nargin<3
   error('Invalid syntax, 3 arguments required');
elseif nargin==3
   TS = 0;
end

if isa(dcval,'double') && ndims(dcval)==2 && ...
      all(size(dcval)==[1 1]) && imag(dcval)==0
   if isa(ginf,'double') && ndims(ginf)==2 && ...
         all(size(ginf)==[1 1]) && imag(ginf)==0
      if (abs(dcval)<1 && 1<abs(ginf)) || ... 
            (abs(dcval)>1 && 1>abs(ginf))
         if isa(bw,'double') && ndims(bw)==2 && ...
               all(size(bw)==[1 1]) && imag(bw)==0 && bw>0 && ~isinf(bw)
           
            if TS==0
               poleloc = abs(bw)*sqrt((ginf*ginf-1)/(1-dcval*dcval));
               G = ss(tf([ginf poleloc*dcval],[1 poleloc]));
            elseif TS>0 || TS==-1
                TSabs = abs(TS);
                if bw*TSabs>pi
                    error('CROSSW*TS must be less than PI.');
                end 

                % Apply bilinear transform to cont.-time solution 
                poleloc = abs(bw)*sqrt((ginf*ginf-1)/(1-dcval*dcval));
                bt = bw*sin(bw*TS)/(cos(bw*TS)-1);
                
                aa = poleloc*dcval-ginf*bt;
                bb = poleloc*dcval+ginf*bt;
                cc = poleloc-bt;
                dd = poleloc+bt;
                G = ss(tf([aa bb],[cc dd],TS));     
                
               %Following is equivalent but uses c2d 
               %poleloc = abs(bw)*sqrt((ginf*ginf-1)/(1-dcval*dcval));
               %Gc = ss(tf([ginf poleloc*dcval],[1 poleloc]));
               %G = c2d(Gc,abs(TS),'prewarp',bw); 
               %if TS==-1
               %    G.Ts = TS;
               %end
            else
               error('Invalid Sampling Time');
            end
            
         else
            error('Invalid CrossOver Frequency');
         end
      else
         error('Invalid relation between the low and high frequency gains.');
      end
   else
      error('Invalid High Frequency gain.');
   end
else
   error('Invalid Low Frequency gain.');
end