findtunetom2.m 2.94 KB
Newer Older
Thomas Perron's avatar
Thomas Perron committed
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29
function [tunef,numdc,peakss,phasef,camp,uvect]=findtunetom2(M,divider,window)

%[tunef,numdc]=findtunetom2(M,divider)
%finds thew tune for real or complex sets of orbits. ALl orbits should have
%the same length. one orbit is one collumn, and there ashould be as many
%collums
%as there are different orbits to treat.
% the size of M should be orbit_lengthxnumber_of _orbits. 
% WINDOW should be equal to 1 if a hanex window should be aplied.

%*******************************************************************
%THE VALUE OF THE TUNE IS GIVEN CONSIDERED THAT THE ORBITS ARE RED ON A
%SINGLE BPM (SAMPELING AT REVOLUTION FREQUENCY). IF THE ORBIT IS READ ON
%MORE THAN ONE BPM THE GIVEN TUNE VALUE SHOULD BE MULTIPLIED BY THE NUMBER
%OF BPMS CONSIDERED.
%EG: ALL 224 BPM ON A FIRST TURN: MULTIPLY THE COMPUTED TUNE
%VALUE BY 224
%*******************************************************************


% tune is then a line vector containing the tunes. The tune is calculated
% doing a hanning window and making an interpolation to find the right peak
% value on the spectrum.
%M is the orbit matrix, DIVIDER is the step to resolve the main peak (Hz). step=
%FFTstep/DIVIDER
phases=NaN;
tune=NaN;
pointnotok=~isfinite(M);
size(M);
Thomas Perron's avatar
Thomas Perron committed
30
indbaddata=sum(pointnotok,1)>size(M,1)*0.8;
Thomas Perron's avatar
Thomas Perron committed
31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108
M(pointnotok)=0;
numdc=0;
c=size(M,2);
n=size(M,1);
nu=[0:1:n-1];
han=(0.53836-0.46164*cos(2*pi*nu./(n-1)))';
hanext=han(:,ones(1,c));
if nargin>2
    if window==1
M=M.*hanext;   
    end
end

res=abs(fft(M(:,~indbaddata)));
Mok=M(:,~indbaddata);
res([1:3 end-3:end],:)=0;
matexp=(i*2*pi*nu./n)';

[bid,ind]=max(res);
if isreal(M) res(round(n/2):n,:)=0; end

for k=1:size(res,2)
    serie=res(:,k);
          
    if ind(k)==1
     subind=[0:1/(divider-1):1];
   col=Mok(:,k);
   
   matcomp=sum(exp(matexp(:,ones(1,divider)).*subind(ones(n,1),:)).*col(:,ones(1,divider)));
    prodscal=abs((matcomp));
    p=phase((matcomp));  
    [bid,indn]=sort(prodscal,2,'descend');
   peakss(k)=prodscal(indn(1)); 
   tune(k)=(subind(indn(1)))/n;
   phases(k)=p(indn(1));
   camp(k)=matcomp(indn(1))/n;
   size(matexp(:,ones(1,divider)))
   size(subind(ones(n,1),:))
   
   uvect(:,k)=exp(matexp.*subind(indn(1)));
   %exp(matexp(:,ones(1,divider).*subind(ones(n,1),:)));
   %peak(k)=0
    %    tune(k)=0;
     %  numdc=numdc+1;
    %elseif ind(k)==n
    %tune(k)=1-1/n;
    
   
    
    else 
    a=ind(k);
    col=Mok(:,k);
    
    %subind=[ind(k)-1.5:1/(divider-1):ind(k)-0.5];
   subind=[ind(k)-1.5:1/(divider-1):ind(k)-0.5];
   
   
   
   matcomp=sum(exp(matexp(:,ones(1,divider)).*subind(ones(n,1),:)).*col(:,ones(1,divider)));
   p=phase((matcomp));  
   prodscal=abs((matcomp));
    [bid,indn]=sort(prodscal,2,'descend');
  
    peakss(k)=prodscal(indn(1)); 
   tune(k)=subind(indn(1))/(n);
    phases(k)=p(indn(1));
     camp(k)=matcomp((indn(1)))/(n);  
   
     uvect(:,k)=exp(matexp.*subind(indn(1)));
    end
  
    
end
    
tunef(indbaddata)=NaN;
   tunef(~indbaddata)=tune;
phasef(indbaddata)=NaN;
phasef(~indbaddata)=phases;