tic
clc
close?all
clear
load?guwu09.txt;
load?guyou09.txt;
load?guyou095.txt;
FONTSIZE=0.3;
TotalNumOfInterations=30;

OriginalSignal=guyou095（51:5622）;
%?OriginalSignal=guyou09（:2）;
plot（OriginalSignal）;?xlabel（‘節拍‘）;?ylabel（‘速度‘）;??hold?on;
set（gca‘xlim‘[1?600]）;
%?OriginalSignal=OriginalSignal‘;
SignalLength=length（OriginalSignal）;%clear?g4;
OriginalSignal（SignalLength+1:2^ceil（log2（SignalLength）））=0;
SignalLength=length（OriginalSignal）;
CurResidues=OriginalSignal;
TotalActualCorre=0;
SamplingFrequency=1;?????????????%the?sampling?frequence?is?application?dependent
Numfrequency=SignalLength;?????%number?of?samples?in?frequency?domain
TotalEnergyOfSignal=sum（CurResidues.*CurResidues）;
Iteration=1;
T=1:SignalLength;
T=T‘;
AverageInterationTime=0;
AverageRefineTime=0;
AtomsCalculateTime=0;
?re=OriginalSignal;
%--------------------------------------------------------------------------
%?main?loop?
%?search?or?the?best?atom?which?match
%?the?main?component?of?the?residues?of?the?signal?uring?each?ineration
%?----------------------------------------------------------------
jmax=round（log2（SignalLength））;%max?number?of?scales?in?the?subset（da）for?searching?0a=2;
DetaU=0.5;
DetaW=pi;
%?TempComGr=zeros（SignalLength1）;
%?TempRealGr=zeros（SignalLengthTotalNumOfInterations）
TempCoef=zeros（TotalNumOfInterations5）;
ATOMS=cell（1jmax-1）;
TOTAL_TIME=cputime;
for?Q=1:jmax-1%max?number?of?translations?in?the?subset?for?searching?0<=p<=Length?2^（-Q+1）
%-------------------------------------------
%?????generate?the?gr?for?current?r（Qpk）
%------------------------------------------
???TimeScale=a^Q
???Timeshift=ceil（SignalLength/2）;
???TimeUp=Timeshift+TimeScale;
???TimeDown=Timeshift-TimeScale+1;
???if?TimeUp>SignalLength?TimeUp=SignalLength;end
???if?TimeDown<=0?TimeDown=1;end
???TimeUpTimeDown（Q:）=[TimeDownTimeUpTimeScale];
???ATOMS{Q}=exp（-pi*（（T（TimeDown:TimeUp）-Timeshift）/（TimeScale））.^2）;
end
AtomsCalculateTime=cputime-TOTAL_TIME;
%?----------------------------------------
TOTAL_TIME=cputime;
while（TotalNumOfInterations>=Iteration）??%main?loop?for?looking?for?best?atom?each?time
????TIME=cputime;
????Residues_FFT=fft（CurResiduesNumfrequency*2）;
????[TempOptimFrequency]=max（（Residues_FFT））;?????%here?we?use?fft?to?select?the?frequency?center?of?the?atoms
????OptimFrequency=（OptimFrequency）/Numfrequency*pi;
?for?Q=1:jmax-1
?????CurGr=ATOMS{Q}.*exp（j*（OptimFrequency）*T（TimeUpTimeDown（Q1）:TimeUpTimeDown（Q2）））;
?????CurGr=CurGr/norm（CurGr）;
?????CurGr（2*TimeUpTimeDown（Q3）+1:SignalLength*2）=0;
?????qweqweq=fft（CurGrSignalLength*2）;
?????CurRffGrcom=ifft（conj（qweqweq）.*Residues_FFT）;
?????CurRffGrcom=[CurRffGrcom（2*SignalLength-TimeUpTimeDown（Q3）+1:2*SignalLength）;CurRffGrcom（1:2*SignalLength-TimeUpTimeDown（Q3）+1）];
?????[CurRffGrPshift]=max（CurRffGrcom）;
?????CurRffGr=abs（CurRffGr）;??????%to?cauculate?the?
?????if?CurRffGr>TempCoef（Iteration1）????%store?