Contents

%generate periodic function with uniform sampling
%add noise
%OUPUPT: y(n) written to output file
%
%Profile: a rising function (integrating capacitor) followed by
%a decaying capacitor discharge

Specify Inputs

    OutFile='PeriodicSignal.dat';
    P=53.7324;  %period in seconds
    dt=0.1;     %sampling interval in seconds
    Nperiod=1000; %number of periods in time series
    f1=0.1;     %fraction of period for the rising part
    mbin=21;    %number of bins for folding
    alpha=1.4;   %want the decay to be faster than exponential

    a=1;    %amplitude of Gaussian noise

    fprintf('\n\n')
    fprintf('Generating exponential charge & discharge profile\n');
    fprintf('Period: %f (s)\n',P);
    fprintf('sampling interval %f (s)\n',dt);
    fprintf('Number of periods of time series %d\n',Nperiod);
    fprintf('Amplitude of Gauassian noise %f\n',a);

Generating exponential charge & discharge profile
Period: 53.732400 (s)
sampling interval 0.100000 (s)
Number of periods of time series 1000
Amplitude of Gauassian noise 1.000000

Generate the basic profile

    t=0:dt:Nperiod*P;
    n=numel(t);

    y=GenerateProfile(t,P,f1,alpha);  %code can be found at bottom of page
    yn=y+a*randn(1,n);

Write out the timeseries

    fid=fopen(OutFile,'w');
    fprintf(fid,'%%time series of a periodic signal\n');
    fprintf(fid,'%% column 1: timme, column 2: amplitude\n');
    fprintf(fid,'%f %f\n',[t;y]);
    fclose(fid);

Fold the time series

    mperiod2plot=4;
    figure(1)
    subplot(2,1,1)
        ind=1:fix(mperiod2plot*P/dt);
        plot(t(ind),y(ind))
        str=sprintf('first %d periods',mperiod2plot);
        title(str)
        ylabel('amplitude')
    subplot(2,1,2)
        plot(t(ind),yn(ind))
        title('with noise')
        xlabel('time(s)')
        ylabel('amplitude')


       dP=0.01;   %fold at period P+dP
       fprintf('true period, P:  %f\n',P);
       fprintf('fold period is P+dP where dP=%f\n',dP);
       fprintf('will fold into %d bins\n',mbin);
       Ptry=P+dP;
       figure(2)
       FoldProfile(t,yn,Ptry,mbin);

   fprintf('\n')

true period, P:  53.732400
fold period is P+dP where dP=0.010000
will fold into 21 bins

GenerateProfile

The function GenerateProfile is given below

function y=GenerateProfile(t,P,f1,alpha)

    t=t-fix(t/P)*P;
    a=f1*P;
    b=P-a;
    ycross=1-exp(-1);
    y=(t<a).*(1-exp(-t/a))+(t>a).*ycross.*exp(-(t-a).^alpha/b);
end

Published with MATLAB® R2016b