clc;
clear all;
close all;


disp('*****************************************');
disp('***        MATLAB - SBPA V0.3         ***');
disp('*****************************************');

%% Param鑤res
vitesse_son = 340;              %340 m/s
distance_robot_balise = sqrt(3^2+2.1^2);      %~3.6 m
f_us = 32e3;
F_osc = 128e6;
p = 1;

%% Longueur et fr閝uence SBPA
t_robot_balise = round(distance_robot_balise / vitesse_son * 1000) / 1000;

t_sbpa_haut = 2.54/f_us;

s = sprintf('Longeur minimum de la SBPA  : %f ms',t_robot_balise*1000);
disp(s);
s = sprintf('Longeur minimum de l''閠at haut de la SBPA : %f ms',t_sbpa_haut*1000);
disp(s);
disp('.');
disp('Syst鑝e d''equations :');
s = sprintf('\t nr x Tsbpa > %f ms',t_sbpa_haut*1000);
disp(s);
s = sprintf('\t (2^nr - 1) x Tsbpa > %f ms',t_robot_balise*1000);
disp(s);

%r閟olution de l'閝uation
temp='';
temp = strcat(temp,num2str(t_sbpa_haut));
temp = strcat(temp,' * 2^(nr) - ');
temp = strcat(temp,num2str(t_robot_balise));
temp = strcat(temp,' * nr = ');
temp = strcat(temp,num2str(t_sbpa_haut));
sol = solve(temp);

nr = round(double(sol(2)));
Tsbpa = t_sbpa_haut / nr;

f_ech = p/Tsbpa;

longueur_sbpa = 2^nr - 1;
nb_ech = round(longueur_sbpa * f_ech * Tsbpa);



%% Calcul des param鑤res du dsPIC

ADCS = 2*(F_osc)/(4*f_ech*13)-1;
ADCS = round(ADCS);
f_ech_pic = 2*F_osc/(4*13*(ADCS+1));

PR1 = F_osc/(4*f_ech_pic);

%% Affichage R閟ultats
disp('*****************************************');
disp('R閟ultats :');
s = sprintf('\t Nombre de registres : \t\t\t    %d registres',nr);
disp(s);
s = sprintf('\t Longeur de la SBPA : \t\t\t  %d 閘閙ents',longueur_sbpa);
disp(s);
s = sprintf('\t Fr閝uence de la SBPA : \t\t% d Hz',1/Tsbpa);
disp(s);
s = sprintf('\t Fr閝uence d''閏hantillonnage : \t %d Hz',f_ech);
disp(s);
s = sprintf('\t Nombre d''閏hantillons n閏essaire pour 1 SBPA : %g 閘閙ents',nb_ech);
disp(s);

disp('*****************************************');
disp('Param鑤res du dsPIC :');
disp('TIMER 1 :');
s = sprintf('\t Valeur de PR1 : \t\t  %d',PR1);
disp(s);
s = sprintf('\t Pr閟calaire du timer 1 : \t%d',1);
disp(s);
s = sprintf('\t Fr閝uence TIMER 1 : \t%d Hz',f_ech_pic);
disp(s);
disp('ADC :');
s = sprintf('\t ADC_SAMPLE_TIME_%d',0);
disp(s);
s = sprintf('\t ADC_CONV_CLK_%dTcy',round(ADCS/2));
disp(s);
disp('*****************************************');
s = sprintf('\t Fr閝uence d''閏hantillonnage du dsPIC : \t\t%g kHz',f_ech_pic*1e-3);
disp(s);
s = sprintf('\t Temps d''閏hantillonnage du dsPIC : \t\t\t  %g 祍',1e6/f_ech_pic);
disp(s);
disp('*****************************************');
disp('.');

%% Calcul charge m閙oire :
mem = 3 * 2 * 2 * nb_ech;

s = sprintf('Charge m閙oire pour trois signaux : %g octets, soit %g ko',mem, mem/1024);
disp(s);

%% Affichage 閏hantillonnage
t_ech = 1/f_ech;

%Chargement des donn閑s
load donnee.mat

%centrage des donn閑s
ech = ech - 512;
%%
nb_ech=1024;
donnee1=ech(1:nb_ech);
donnee2=ech(nb_ech+1:2*nb_ech);
donnee3(1:100)=donnee1(nb_ech-99:nb_ech);    %test de corr閘ation
donnee3(101:nb_ech)=donnee1(1:nb_ech-100);

correlation=xcorr(donnee1, donnee2);
x_corr = [-nb_ech+1:nb_ech-1];
figure('Name', 'Convolution', 'NumberTitle','off');
subplot(2,1,1);
plot(x_corr,correlation);
title('Corr閘ation de donnee1 et donnee2');

correlation2=xcorr(donnee1, donnee3);
x_corr = [-nb_ech+1:nb_ech-1];
%figure('Name', 'Convolution2', 'NumberTitle','off');
subplot(2,1,2);
plot(x_corr,correlation2);
title('Corr閘ation de donnee1 et donnee3');

%coupure de ech
if length(ech) > nb_ech
    ech = ech(1:nb_ech);
end
%%
% donnee1=ech(1:nb_ech);
% donnee2=ech(nb_ech+1:2*nb_ech);
% correlation=xcorr(donnee1, donnee2);
% x_corr = [-nb_ech+1:nb_ech-1];
% figure('Name', 'Convolution', 'NumberTitle','off');
% plot(x_corr,correlation);
% 
% %coupure de ech
% if length(ech) > 2*nb_ech
%     ech = ech(1:2*nb_ech);
% end
%%
disp('.');
s = sprintf('Taille de l''閏hantillonnage : \t\t%g 閘閙ents',length(ech));
disp(s);
s = sprintf('Valeur moyenne de l''閏hantillonnage : \t\t%f',mean(ech));
disp(s);
disp('.');

%Vecteur temps
x = [0:t_ech:(size(ech)-1)*t_ech];

%affichage
figure('Name', 'Echantillonnage', 'NumberTitle','off');
plot(x,ech);
title('Echantillonnage', 'fontsize', 16);
xlabel('Temps en s');

%affichage comparatif
figure('Name', 'Echantillonnage comparatif', 'NumberTitle','off');
subplot(3,1,1);
plot(donnee1);
title('Donnee1');
subplot(3,1,2);
plot(donnee2);
title('Donnee2');
subplot(3,1,3);
plot(donnee3);
title('Donnee3');

%affichage comparatif
figure('Name', 'Echantillonnage comparatif', 'NumberTitle','off');
hold on;
plot(donnee1,'color',[1 0 0]);
plot(donnee2,'color',[0 0 1]);


%% Cr閍tion de la SBPA
reg_dec = [1,0,1,0,1,0,1,0,1,0,1]; %SBPA de d閜art
freq = 136e3;
k=0:1/freq:(longueur_sbpa-1)/freq;
f = freq*linspace(0,1,longueur_sbpa);
%k=1:longueur_sbpa;        %Vecteur pour affichage SBPA

for C = 1:2*longueur_sbpa
    
    %les ou exclusive
    temp=xor(reg_dec(1,9),reg_dec(1,11));
    
    %le decalage
    reg_dec(:,[2:11]) = reg_dec(:,[1:10]);
    reg_dec(1,1) = temp;
    
    SBPA(1,C) = reg_dec(1,1);
    
end

%centre en 0 la SBPA
for C = 1:2*longueur_sbpa
   if SBPA(1,C) == 0;
       SBPA(1,C)=-1;
   end
end

%****************************
% AFFICHAGE DE LA SBPA
%****************************
SBPA=SBPA(1,[1:longueur_sbpa]);
figure('Name','SBPA', 'NumberTitle','off');
stairs(k,SBPA);

spectre = abs(fft(SBPA));
figure('Name','Spectre SBPA', 'NumberTitle','off');
plot(f,spectre);

ech3=ech3([1:longueur_sbpa]);
ech3 = ech3 - mean(ech3);

spectre = fft(ech3);
figure('Name','Spectre ech', 'NumberTitle','off');
plot(f,spectre);