Discussion :

[d0n32]

BOnjour tous,
Nouvelle journée, nouveau problème: je souhaite mettre en place l'algorithme de calcul d'une fonction trigo, via le fameux algorithme CORDIC (celui implémenté dans les calculatrices).
Celui-ci ne fonctionne que pour une valeur d'angle comprise entre 0 et pi/2
Voici ce que j'ai écrit :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

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 v = cordic(beta,n)
 
%Initialisation
v=[1;0];
sigma=1;
Kn = prod(1./sqrt(1+2.^(-2*(0:(n-1)))));
 
if (beta>=pi && beta<=3*pi/2)
    cordic(beta-pi,n);
    v=-v;
    return;
end
 
%Algorithme Cordic pour beta entre 0 et pi/2 (algorithme non valable sinon)
if (beta>=0 && beta<=pi/2)
    if beta==pi/4  % Bug dans algo Cordic sinon
        v=[sqrt(2)/2;sqrt(2)/2];
        return;
    end
    %Itérations au calcul de cos(beta) et sin(beta)
    for i=0:n-1
        R=[1 -sigma*(2^(-i)); sigma*(2^(-i)) 1];
        v=R*v;
        beta = beta-sigma*atan(2^(-i));
        sigma = sign(beta);
    end
    %Calcul final du vecteur (cos(beta),sin(beta))
    v=v*Kn;
end

Je souhaiterais étendre l'algorithme à des angles compris en tre 0 et 2*pi. Je commence à traiter le cas (pi, 3pi/2), cf mon code.
Le problème est que quand je teste un truc du style cordic(5*pi/4), j'obtiens un truc incohérent (juste -v, comme si le rappel de la fonction ne marchait pas).
Une idée SVP ?
Merci d'avance

[Invité]

Bonjour,

C'est normal, tu ne donnes pas le résultat de cet appel à v

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
if (beta>=pi && beta<=3*pi/2)
    v =  - cordic(beta-pi,n);
    return;
end

[d0n32]

je te remercie encore ! J'ai un nouveau problème lié à l'affectation de variable :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
if (beta>pi/2 && beta<pi)
    v=cordic(beta-pi/2);
    v(1,:)=-v(2,:);
    v(2,:)=v(1,:);
    return;
end

J'aimerais que pour v(2,: )=v(1,: ); ce soit l'ancien v(1,: ) qui soit pris en compte, ce qui n'est pas possible dans cette configuration. Mais je ne vois pas comment m'en dépatouiller..


EDIT :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
if (beta>pi/2 && beta<pi)
    v=cordic(beta-pi/2);
    a=v(1,:);
    v(1,:)=-v(2,:);
    v(2,:)=a;
    return;
end

contourne le problème !

[Jean Dumoncel]

Bonjour,

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
    a=v(1,:);
    v(1,:)=-v(2,:);
    v(2,:)=a;

s'écrit aussi :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

v(:,1:2) = [v(:,2) v(:,1)];

ou encore :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

v(:,1:2) = v(:,2:-1:1);

[d0n32]

Merci

[d0n32]

Bonjour !
Désolé de déterrer le sujet, mais en me repenchant sur mon algorithme, je me suis aperçu de son extrême lenteur ! Devant l'appliquer plusieurs fois de suite à un gros paquets de données (125k pour être précis), Matlab met plus de 4 minutes à trouver la réponse !
J'ai diminué le nombre d'itérations à 15 pour limiter le nombre de boucles. JH'ai fait plusieurs tests avec un nombre d'itérations différent, et le résultat est surprenant : avec seulement 3 itérations, il met juste 2 secondes de moins pour calculer, qu'avec 15 itérations ! J'en conclus donc que la lenteur ne vient pas de la boucle elle-même..?
J'ai beau chercher, je ne vois pas.. Auriez-vous une idée ? Si besoin je peux reposter le code, mais il est déjà présent un peu plus haut dans le topic.. La seule différence est que j'ai posé n=15, et non plus une variable de la fonction.
Merci d'avance

[d0n32]

Bon je n'ai pas trouvé la réponse mais à ma question, mais j'ai essayé de remplacer mon algorithme par une implémentation direct du cosinus/sinus. J'ai écrit un truc du genre

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
 
for k=1:length(y)
 x(1,k)=cos(y(k));
 x(2,k)=sin(y(k));
end

qui est équivalent à l'utilisation de l'agorithme Cordic, et Matlab a à peine mis 3 secondes de moins à trouver le résultat.
J'en déduis que c'est juste lié au nombre de valeurs que je demande à calculer, et qu'à moins de diminuer le nombre de points, je n'arriverai pas à gagner beaucoup de temps en calcul !

[Jean Dumoncel]

Bonjour,

peux-tu nous remontrer le code complet que tu executes? et un exemple de l'appel de cette fonction?

[d0n32]

Bonjour,
Voici la fonction :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

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
30
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
function [v] = cordic(beta)
%Initialisation
n=15;  %Nombre d'itérations
v=[1;0];
sigma=1;
Kn = prod(1./sqrt(1+2.^(-2*(0:(n-1)))));
a=0; b=0; %Variables intermédiaires
 
%Algorithme Cordic pour beta entre 0 et pi/2 (algorithme non valable sinon)
if (beta>=0 && beta<=pi/2)
    if beta==pi/4  % Bug dans algo Cordic sinon
        v=[sqrt(2)/2;sqrt(2)/2];
        return;
    end
    %Itérations au calcul de cos(beta) et sin(beta)
    for i=0:n-1
        R=[1 -sigma*(2^(-i)); sigma*(2^(-i)) 1];
        v=R*v;
        beta = beta-sigma*atan(2^(-i));
        sigma = sign(beta);
    end
    %Calcul final du vecteur (cos(beta),sin(beta))
    v=v*Kn;
    return;
end
 
 
%pour beta entre pi/2 et pi
if (beta>pi/2 && beta<pi)
    v=cordic(beta-pi/2);
    a=v(1,:);
    v(1,:)=-v(2,:);
    v(2,:)=a;
    return;
end
 
%pour beta entre pi et 3*pi/2
if (beta>=pi && beta<=3*pi/2)
    v = -cordic(beta-pi);
    return;
end
 
%pour beta entre 3*pi/2 et 2*pi
if (beta>3*pi/2 && beta<=2*pi)
    v=cordic(beta-3*pi/2);
    b=v(1,:);
    v(1,:)=v(2,:);
    v(2,:)=-b;
    return;
end
 
if beta==-pi/2  % Seul cas où phi est négatif (cf phi_br=+/- pi/2)
    v=[0;-1];
    return;
end

et un exemple d'application de la fonction :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
for k=1:length(phase)
    y2(:,k)=cordic(phase(k));
end

où phase est un signal de 125000 points

[d0n32]

Bon, je me suis résigné à devoir attendre que Matlab calcule si lentement.
J'ai par contre un autre problème, que je trouve assez étrange :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
 
%Initialisation
fe=125e6;
t_simu=10^-3;
fs=70e6;
t=linspace(0,1e-3,fe*t_simu);
 
%Signal d'input
y=[cos(2*pi*fs*t+pi/4);
    sin(2*pi*fs*t+pi/4)];
 
%Signaux de sortie du cordic
phase=nco(fe,fs,t_simu);
for k=1:length(phase)
    y2(:,k)=cordic(phase(k));
end
%% Multiplication des signaux y et y2
s_avCIC(1,:)=2*y(1,:).*y2(1,:);
s_avCIC(2,:)=2*y(2,:).*y2(1,:);

A priori, nco+cordic fournit un signal sinusoïdal. Et cos(a)*cos(b) = 0.5*(cos(a+b) + cos(a-b)) donc dans mon cas, je suis censé obtenir en sortie un signal sinusoïdal de fréquence 2*fs, centré autour de la valeur cos(pi/4) (ou sin(pi/4) pour la deuxième coordonnée).
Le problème est que ce n'est pas du tout ce que j'obtiens !
Je devrais avoir ça : une sinusoïde de HF centrée autour de cos(pi/4) = 0.707


J'obtiens ça : une sinusoïde de HF qui oscille autour d'un signal BF (manifestement)..


J'ai checké le cordic (cf code sur le post précédent), et le souci ne vient pas de lui. Je ne vois donc que le NCO !
Voici son code, largement inspiré des explications fournies sur Wikipédia :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
 
function [phase] = nco(fe,fi_br,t_simu)
 
phase = zeros(1,t_simu*fe); %Phase
increment = 2^48*fi_br/fe;
 
for j=1:length(phase)-1
    phase(j+1)=mod(phase(j)+increment,2^48-1);
end
phase=phase*2*pi/2^48;

[d0n32]

Le problème pourrait-il venir de linspace ? En entrant

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

t=linspace(0,1e-3-8e-9,t_simu*fe);

au lieu de

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

t=linspace(0,1e-3,t_simu*fe);

, les résultats sont ceux escomptés !