% ion1 is an example of an ion distribution function.

ALLpole =[];
ALLresid=[];


Trat=100;			% Trat=Te/Ti Temperature ratio

n_vec    =[5   0];		% # of terms for Maxwell expansion
nb_vec   =[0   0];		% # of terms for Butterworth 
vtbut_vec=[2.0 0];		% Cut velocity in cs units
op2      =[1   0];		% densities
vt       =[1   1]/sqrt(Trat);	% thermal velocity cs==1
vd       =[0   0];		% drift velocity cs==1

op2=op2/sum(op2);		%use total plasma frequency
%%%% Normalize  
vtbut_vec=vtbut_vec./vt;
for jj=1:1
  n= n_vec(jj);  
  nb=nb_vec(jj);
  vtbut=vtbut_vec(jj);

  pol=0*(1:2*n+1);
  pol(2*n+1)=1;
  for i=1:n				% Poles for Maxwellian
    pol(2*(n-i)+1)=1.0/(2^i*prod(1:i));
  end
  sol=roots(pol);			% These are my poles
  res=pole2residue(sol);		% Find the residues
 
  z0=butterworthpoles(nb,vtbut); 	% Butterworth poles


  allpole =[sol;z0];			%poles from expansion and mask
  allresid=pole2residue(allpole);	% Find residues
  [uhppole,uhpresid]=uhp(allpole,allresid); % Keep poles in upper half plane

  allresid=allresid/(2i*pi*sum(uhpresid)); % normalize zeroth moment
  for j=1:length(allpole)                  % return to variables with dimension
    ALLpole =[ALLpole allpole(j)*vt(jj)+vd(jj)];
    ALLresid=[ALLresid op2(jj)*allresid(j)];
  end
  
end

[UHPpole,UHPresid]=uhp(ALLpole,ALLresid); % Keep poles in upper half plane

