Algorithme CORDIC, souci..

Version imprimable

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:

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 :ccool::mrgreen:

17/06/2010, 12h44
Invité
Bonjour,

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

1 2 3 4 if (beta>=pi && beta<=3*pi/2) v = - cordic(beta-pi,n); return; end
17/06/2010, 13h28
d0n32
je te remercie encore ! J'ai un nouveau problème lié à l'affectation de variable :
Code:

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..
:ccool:

EDIT :
Code:

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 !
17/06/2010, 13h50
Jean Dumoncel
Bonjour,
Code:

1 2 3 a=v(1,:); v(1,:)=-v(2,:); v(2,:)=a;
s'écrit aussi :

Code:

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

ou encore :

Code:

v(:,1:2) = v(:,2:-1:1);
17/06/2010, 14h48
d0n32

Merci :mrgreen:
07/07/2010, 09h34
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 :ccool:
07/07/2010, 14h01
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:

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 !
07/07/2010, 14h17
Jean Dumoncel

Bonjour,

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

Bonjour,
Voici la fonction :

Code:

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:

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

où phase est un signal de 125000 points

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:

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
http://img188.imageshack.us/img188/9094/12186490.jpg

J'obtiens ça : une sinusoïde de HF qui oscille autour d'un signal BF (manifestement)..
http://img688.imageshack.us/img688/9779/69600207.jpg

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:

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;

09/07/2010, 13h19
d0n32

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

Code:

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

au lieu de

Code:

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

, les résultats sont ceux escomptés !
09/07/2010, 17h13
d0n32

Bon, en fait j'ai toujours le problème, dès que j'essaie de répéter plusieurs fois de suite le code donné 2 posts plus haut en changeant cette ligne

Code:

phase=nco(fe,fs,t_simu);

par

Code:

phase=nco(fe,fs+k,t_simu);

où k est incrémenté à chaque fois que je répète l'algorithme, je me retrouve au final avec des courbes similaires à celle affichées dans le précédent post.. C'est assez énervant..
Si vous avez une idée, of course, je suis preneur :ccool:
12/07/2010, 09h35
d0n32

Autre phénomène étrange : le premier tour de boucle du COrdic est réalisé en 95 secondes environ, et les suivants en à peu près 12 secondes, soit 8 fois moins de temps ! Cela est-il normal ?
13/07/2010, 16h48
d0n32

Pour répondre à moi même et clore le sujet, j'ai testé plusieurs types de boucles, en les répétant systématiquement, et j'ai pu observer que Matlab met bien plus de temps à calculer la boucle la première fois, que les fois d'après !
J'en déduis donc que c'est normal.
14/07/2010, 00h13
Invité

Bonsoir,

Je te conseille de regarder la :faq:
15/07/2010, 09h21
d0n32

Salut,
Tiens c'est marrant, je n'avais pas vu les choses comme ça. Cependant, j'avais lu le message de Matlab "considerate preallocating for speed" ou quelque chose du genre, et j'avais initilalisé, comme il en est question dans la FAQ. Et en chronométrant, il n'y avait pas de différences notables ! J'ai donc supprimé la préallocation.. Je vais la laisser, on ne sait jamais :D