Spectre sous Matlab

Résolu/Fermé
Pep7 - 11 mars 2009 à 15:28
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 - 13 mars 2009 à 18:17
Bonjour,
j'essaye de représenter le spectre d'un signal d'après sa formule, mais je n'arrive pas.
Pouvez vous m'aider.
Voici la formule:

l_b (ω)=Q_b∙exp(-(ω^2∙σ_b^2)/2)∙∑_(k=1)^(N_b)▒e^(-jωkT_b )

Avec Q_b = 2.33 * 10^ -9
σ_b = 0.5 *10^ -3
N_b = 210
T_b = 666 * 10^ -12

Merci de votre réponse
Voilà ou j'en suis, mais ça ne joue pas

% Temps
fe = 3*10^6;
Te = 1/fe;
N = 1000;
t = (0:N-1)*Te;
w = 2*pi./t;

Qb = 2.33*10^-9;
Nb = 210;
ob = 0.5*10^-3;
Tb = 666*10^-12
Somme = 0
for u=1:N
for k=1:Nb
Somme = Somme + exp(w(u)*k*Tb)
end
Ib(u) = Qb*exp(-((w(u)^2*ob^2)/2))*Somme
end

figure(8)
plot(w,Ib);grid on,
A voir également:

20 réponses

Char Snipeur Messages postés 9696 Date d'inscription vendredi 23 avril 2004 Statut Contributeur Dernière intervention 3 octobre 2023 1 297
12 mars 2009 à 12:10
Sacabouffe a raison, tu t'y prend mal avec les gens qui te donne des conseils.
Ta somme ne sert pas à rien, il est inutile de faire une boucle pour la calculer car c'est une suite géométrique et qu'il existe une relation pour la calculer simplement (cf. cours du lycée).
Par contre, je ne comprends pas, dans ton équation, tu met "-j" qui est je suppose le nombre imaginaire (j²=-1) mais tu ne le met pas dans ton équation sous matlab.
Sinon, tu pourrai donner plus de détails que "ça ne joue pas" (d'ailleurs, dans le genre anglicisme à la con, c'est pas mal), je suis sur que si tu as mal écris quelque chose, MATLAB râle plus que ça.
5
Merci char sniper.
Désolé pour avant.

oui, j'ai oublié de mettre le j.

Comment la calculer la suite géométrique si j'ai deux paramètres avec différentes tailles (k et w).
Car je dois faire la somme avec toutes les valeurs de k pour une seule valeur de w et répéter l'opération pour chaque valeurs de w.

Ce que je dois obtenir, c'est un spectre avec une raie en DC qui rebondi puis une répétition de cette raie à 15Ghz puis 30ghz ... avec chaque fois une diminution de l'amplitude proportionnelle à exp(-w^2*ob^2)/2

Voilà la fonction

wmin = 1
wmax = 50e9
Npts = 1000
dw = (wmax - wmin)/Npts
w = wmin :dw :wmax-dw ;

Qb = 2.33*10^-9;
Nb = 210;
ob = 0.5e-3;
Tb = 666e-12
k=(1:210)

Ib = Qb*exp(-((w.^2*ob^2)/2))*sum(exp(-i*k*Tb*w))

figure(8)
plot(w,Ib);grid on,

merci pour ta réponse

A bientôt
0
Char Snipeur Messages postés 9696 Date d'inscription vendredi 23 avril 2004 Statut Contributeur Dernière intervention 3 octobre 2023 1 297 > Pep7
12 mars 2009 à 15:31
cf. ce qu'a mis sacabouffe. Il a même intégrer le "j" que tu avais oublié je pense (d'où les sinus, car les sommes d'exponentielles imaginaire rendent souvent des fonctions trigo).
0
Pep7 > Char Snipeur Messages postés 9696 Date d'inscription vendredi 23 avril 2004 Statut Contributeur Dernière intervention 3 octobre 2023
13 mars 2009 à 11:32
Hello,

ok pour la somme, je la retrouve.
par contre je ne comprends pas pourquoi tu divise par w, dans la relation de la somme, je trouve la même chose, mais en multipliant par w

Merci
0
Char Snipeur Messages postés 9696 Date d'inscription vendredi 23 avril 2004 Statut Contributeur Dernière intervention 3 octobre 2023 1 297
13 mars 2009 à 15:00
Je ne comprends rien à tes histoires de cos ou sin, mais je trouve comme sacabouffe pour la somme ( https://fr.wikipedia.org/wiki/Suite_g%C3%A9om%C3%A9trique )
sum(k=1,Nb,exp(-ikwT))=(exp(-iwT)-exp(-iwT(Nb+1)))/(1-exp(-iwT))
La partie du dessous donne : exp(-iwT/2)*-2i.sin(wT/2)
La partie du dessus donne : exp(-iwT)(1-exp(-iwTNb))=exp(-iwT)*exp(-iwTNb/2)(2i.sin(wTNb/2))
Soit au final :
-exp(-iwT(NB/2+1/2))·sin(wTNb/2)/sin(wT/2)
Voilà.
2
Char Snipeur Messages postés 9696 Date d'inscription vendredi 23 avril 2004 Statut Contributeur Dernière intervention 3 octobre 2023 1 297
13 mars 2009 à 15:05
Pourquoi poses tu des questions si c'est pour ne pas écouter la réponse ?
As tu au moins testé sur ton PC la solution donnée au premier message ?
Je comprends que tu veuilles comprendre ce qui fait, mais tu n'emploi pas la bonne méthode.
Prend en considération ce que te dises les gens (surtout si tu leur demandes conseil), accepte que tu puisse te tromper.
2
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
13 mars 2009 à 15:29
On a tiré le gros lot là...
Ben ouais mais bon... quand quelqu'un utilise un logiciel de calcul numérique à 500$ sans les toolbox, je pars du principe qu'il connaît le programme de maths de collège et de lycée, je me suis trompé et je te présente mes plus humbles excuses.
Cela dit, ta question était comment implémenter la formule que t'as écrite parce qu'il y avait une coquille dans ton code et la réponse est au message 1 depuis hier ; le premier code correspond à ce que t'as écrit où le calcul de la somme de la série géométrique est fait par Matlab, le second code correspond au calcul de la série géométrique fait la main. Et puis bon, on est dans la forum Programmation, pas dans la cours de récré à attendre la fin du cours de maths pour aller manger à la cantine.
Sans compter que si t'avais montré un minimum de respect pour les personnes qui viennent t'aider au lieu de balancer un truc du genre casse-toi guignol, y a pas quelqu'un d'autre qui veut me répondre ?
bien sur qu'elle sert ma somme.

quelqu'un d'autre peu m'aider
Ben j'aurais sûrement été bien plus patient. D'autant plus que si effectivement, ta somme servait tant que ça, il te suffisait d'utiliser le premier code et pas le deuxième, c'est tout. D'une manière ou d'une autre, le problème était résolu, point.

Repartons donc sur des bases saines et recommençons depuis le début, voudrais-tu tout d'abord, je te prie, aller te documenter sur ce site :
https://www.mathematiquesfaciles.com/
Pour qu'on puisse résoudre ensuite ton problème.

Bon et sinon, ouaip, il manque le 1/2 dans la première exponentielle, alors je vais te montrer comment on le change, tu remplaces exp(-ob^2*w.^2) par exp(-ob^2*w.^2/2), quel que soit le code que tu choisis, voilà...
2

Vous n’avez pas trouvé la réponse que vous recherchez ?

Posez votre question
Char Snipeur Messages postés 9696 Date d'inscription vendredi 23 avril 2004 Statut Contributeur Dernière intervention 3 octobre 2023 1 297
13 mars 2009 à 15:57
J'ai copié le programme de Sacabouffe dans Octave en corrigeant le /2 manquant dans la première exponentielle.
ça tourne très bien. J'obtiens une courbe en cloche, une espèce de demi gaussienne.
Après, ce n'est pas ça que tu voulais obtenir ??? Ba c'est ce que donne la fonction que tu nous a soumis. Es tu sur qu'elle est bonne ta fonction ?
Quant au test, tu met :
Ce que je dois obtenir, c'est un spectre avec une raie en DC qui rebondi puis une répétition de cette raie à 15Ghz puis 30ghz ... avec chaque fois une diminution de l'amplitude proportionnelle à exp(-w^2*ob^2)/2 
Je suppose que c'est de ça dont tu parles lorsque tu dis que tu as testé le code. Mais franchement, encore une fois c'est pas clair.
Comme "ça ne joue pas" qu'est-ce que ça veux dire ??? Il faut bien distingué le problème de langage (error machin...) de l'erreur de logique (tu n'obtiens pas la courbe attendu).
2
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
13 mars 2009 à 16:19
Nan mais de toute façon, c'est n'importe quoi cette formule, on peut pas avoir de quelconque raie à 15GHz ou 30GHz, faut arrêter de fumer la moquette...
La somme est bornée par Nb, donc par 210, quant à Qb*exp(-(w^2*ob^2)/2) pour une fréquence angulaire de 2*pi*15*10^9 donc une fréquence de 15GHz, son ordre de grandeur c'est 1/(10^(4*10^14)), ça fait belle lurette qu'on est passé en-dessous de la précision machine.
Déjà à 15kHz, on est passé en-dessous de la précision machine...
Alors au lieu de brailler comme un âne, il ferait mieux de vérifier sa formule.
2
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
11 mars 2009 à 16:39
Salut
fe = 3*10^6;
w = 2*pi*fe*(0:10^-9:0.5*10^-3);
Qb = 2.33*10^-9;
Nb = 210;
ob = 0.5*10^-3;
Tb = 666*10^-12;
k=(1:Nb).';
Ib=Qb*exp(-ob^2*w.^2).*sum(exp(repmat(Tb*w,size(k)).*repmat(-1i*k,size(w))));
subplot(3,1,1);plot(w,real(Ib));title('Real part');
subplot(3,1,2);plot(w,imag(Ib));title('Imaginary part');
subplot(3,1,3);plot(w,abs(Ib));title('Modulus');
Mais ta somme sert strictement à rien, c'est juste la somme des termes d'une suite géométrique. Donc ce truc suffit amplement.
fe = 3*10^6;
w = 2*pi*fe*(0:10^-9:0.5*10^-3);
Qb = 2.33*10^-9;
Nb = 210;
ob = 0.5*10^-3;
Tb = 666*10^-12;
k=(1:Nb).';
Ib=Qb*exp(-ob^2*w.^2).*exp(-1i*(Nb+1)/2*Tb*w).*sin(Nb*Tb/2*w)./sin(Tb/2*w);
subplot(3,1,1);plot(w,real(Ib));title('Real part');
subplot(3,1,2);plot(w,imag(Ib));title('Imaginary part');
subplot(3,1,3);plot(w,abs(Ib));title('Modulus');
avec un souci quand w*Tb/2≡0[2π], mais bon... suffit de prendre un w pour que ça arrive pas, sinon la somme fait Nb, c'est tout.

Ciao
1
Comprends pas, c'est fait à l'arrache.
bien sur qu'elle sert ma somme.

quelqu'un d'autre peu m'aider

merci
0
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
12 mars 2009 à 11:53
quelqu'un d'autre peu m'aider
Je sais même pas pourquoi je me fais suer à répondre et à donner des lignes de code...
T'as raison, attends quelqu'un d'autre...
Allez, de rien
Ciao
1
Char Snipeur Messages postés 9696 Date d'inscription vendredi 23 avril 2004 Statut Contributeur Dernière intervention 3 octobre 2023 1 297
13 mars 2009 à 12:55
soit plus clair !
1
Voila la forme de sacabouffe

Ib=Qb*exp(-(ob^2*w.^2)/2).*exp(-1i*(Nb+1)/2*Tb*w).*sin(Nb*Tb/2*w)./sin(Tb/2*w)

et moi j'obtiens celle la

Ib=Qb*exp(-(ob^2*w.^2)/2).*exp(-1i*(Nb+1)*Tb*w/2).*sin(Nb*Tb*w/2)./sin(Tb*w/2)

Il y a que le w qui change.
0
Pep7 > Pep7
13 mars 2009 à 13:31
Voilà ce que j obtiens en faite:

Ib=Qb*exp(-(ob^2*w.^2)/2).*exp(-1i*Nb*Tb*w/2).*sin((Nb+1)*Tb*w/2)./sin(Tb*w/2)
0
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
13 mars 2009 à 13:29
Priorité des opérations, programme de collège de 5ème, pas accessible à tout le monde, faut bien le reconnaître...
a/b*c=(a/b)*c=(a*c)/b≠a/(b*c)
1
ok parfait merci par contre encore une petite erreur pour le retour à l'exp, c'est le terme en Nb*Tb*w/2 qui est en cos et en sin.

C'est juste?

merci de ton aide
0
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
13 mars 2009 à 13:38
Non, c'est faux.
Ça fait 24h que le problème est résolu dès le message 1, y a le code avec et sans calcul de la somme, ils donnent les mêmes résultats, alors palabrer inutilement pendant encore 15 jours, ça m'intéresse pas.
Ciao
1
Je suis pas tout a fait d'accord

Tu me propose

Ib=Qb*exp(-(ob^2*w.^2)/2).*exp(-1i*(Nb+1)/2*Tb*w).*sin(Nb*Tb/2*w)./sin(Tb/2*w)

et moi je dis cela

Ib=Qb*exp(-(ob^2*w.^2)/2).*exp(-1i*Nb*Tb*w/2).*sin((Nb+1)*Tb­*w/2)./sin(Tb*w/2)

le terme en sin((n+1)Tb w/2) et le terme en sin(w/2) sont que en sinus
alors que le terme en sin(n Tb w /2) est en cos et en sin selon la somme, donc quand on revient à l'exponentielle,
c'est le terme en cos et sin, soit sin(n Tb w /2)

NON???

Merci
0
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
13 mars 2009 à 14:28
Ouaip...
1
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
13 mars 2009 à 14:33
Ben puisque la somme part de 1, la formule écrite au message 14 est fausse.
On peut y passer encore tout l'après-midi ou tout simplement utiliser les deux lignes de code du message 1...
Ib_by_sum=Qb*exp(-ob^2*w.^2).*...
sum(exp(repmat(Tb*w,size(k)).*repmat(­-1i*k,size(w))));
Ib_direct=Qb*exp(-ob^2*w.^2).*...
exp(-1i*(Nb+1)/2*Tb*w).*sin(Nb*Tb/2*w­)./sin(Tb/2*w);
difference=max(max(abs(Ib_by_sum-Ib_direct)));
>> difference

difference =

  2.1176e-021
1
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
13 mars 2009 à 15:12
Ohlalalala...
Ça dépasse tout ce que j'ai pu imaginer
On a repoussé les limites là
La classe mondiale
Peut-être même le champion du monde

Bon, allons-y avec le parpaing... Maple 12, attention les mirettes...
> simplify(sum(exp(-I*k*Tb*w), k = 1 .. Nb));

                              exp(-I Tb w Nb) - 1
                            - ------------------
                               -1 + exp(I Tb w)
1
Sacabouffe Messages postés 9427 Date d'inscription dimanche 19 août 2007 Statut Membre Dernière intervention 29 mai 2009 1 832
13 mars 2009 à 18:17
Pour terminer cet épisode malheureux de déchéance intellectuelle et faisant abstraction quelques instants du fait qu'à cause de la précision machine, y aura aucune raie, j'ajouterai que la somme géométrique qui a été écrite peut pas donner des raies à 15GHz, 30GHz... Les raies données par une somme de ce type sont situées aux multiples de la fréquence 1/Tb (c'est-à-dire de la fréquence angulaire 2*pi/Tb), c'est à dire ici, aux multiples d'environ 1.5GHz.
Donc 1.5GHz, 3GHz...

Comme l'ordre de grandeur de l'exponentielle est ridiculement petit (en dessous de la précision machine) à la première raie (1.5GHz), et vu que notre ami a pas l'air très futé, il est pas invraisemblable qu'il ait fait une erreur dans les unités. Comme il parle de GHz, ben je dirais qu'il a laissé σ_b exprimé en ns (nanosecondes) ; donc exprimé en secondes, faut le multiplier par 1e-9.

Pour une fréquence d'échantillonnage de 1MHz sur la plage de fréquence [0GHz,15GHz] (10 raies), on a ça...
fe = 1e6;
w = 2*pi*fe*(0:1:15e3);
Qb = 2.33e-9;
Nb = 210;
ob = 0.5e-12;
Tb = 666e-12;
k=(1:Nb).';
Ib=Qb*exp(-ob^2*w.^2/2).*exp(-1i*(Nb+1)/2*Tb*w).*...
sin(Nb*Tb/2*w)./sin(Tb/2*w);
figure(1);plot(w/(2*pi),real(Ib));title('Real part');
figure(2);plot(w/(2*pi),imag(Ib));title('Imaginary part');
figure(3);plot(w/(2*pi),abs(Ib));title('Modulus');
10 raies spectrales : partie réelle
10 raies spectrales : partie imaginaire
10 raies spectrales : module

Pour voir la décroissance, faut utiliser une plage plus large comme [0GHz,1.5THz] (1000 raies), avec une fréquence d'échantillonnage de 1.5MHz, ça donne ça...
fe = 1e6;
w = 2*pi*fe*(0:1.5:1.5e6);
Qb = 2.33e-9;
Nb = 210;
ob = 0.5e-12;
Tb = 666e-12;
k=(1:Nb).';
Ib=Qb*exp(-ob^2*w.^2/2).*exp(-1i*(Nb+1)/2*Tb*w).*...
sin(Nb*Tb/2*w)./sin(Tb/2*w);
figure(1);plot(w/(2*pi),real(Ib));title('Real part');
figure(2);plot(w/(2*pi),imag(Ib));title('Imaginary part');
figure(3);plot(w/(2*pi),abs(Ib));title('Modulus');
1000 raies spectrales : partie réelle
1000 raies spectrales : partie imaginaire
1000 raies spectrales : module

C'est cool le forum Programmation, les gens viennent poser leur torchon incompréhensible et complètement faux, faut deviner la question, faut la résoudre et on se fait insulter comme un malpropre parce qu'ils ont toujours raison et qu'ils sont les plus intelligents du monde... pourquoi ils viennent poser leur question ici alors ? Mystères et boules de gomme...
1
Est ce juste?

∑_(k=0)^(n)▒sin(kx)= sin((n+1)x/2) / sin(x/2) * sin(nx/2)
0
donc ce que je dis au message 14 est juste?

par contre autre problème, cette formule est pour une somme partant de k00, moi elle part de k=1, je dois juste soustraire 1 à la somme totale?

merci
0
vu que e^0 = 1
0
Si la formule 15 est juste, ta formule finale est fausse, mais bon.
En plus tu oublie le /2 de la 1ere expo.

Laisse tombé.

merci quand même d'avoir essayé

Et si jamais la différence est trop importantes pour la mesure de position (au um) d'un faisceau dans des fréquences du GHz
0
Char Snipeur Messages postés 9696 Date d'inscription vendredi 23 avril 2004 Statut Contributeur Dernière intervention 3 octobre 2023 1 297
13 mars 2009 à 15:16
Il faut voir que tu fait du numérique, les erreurs d'arrondi sont obligatoire.
Mais en tout cas, une chose est sur pour diminuer les erreurs, il faut faire un maximum de mathématique, c'est à dire réduire le nombre d'opération réalisées par l'ordinateur. Sachant ça, il vaut mieux calculer la somme à la main comme le fait Sacabouffe que de faire une boucle qui va faire un nombre important d'opération et au final cumulé les erreurs d'arrondis.
Cette remarque est valable pour tout un tas d'opération dont les plus marquante sont intégration et dérivation.
0
Oui j'ai testé, cf lire message 5, avec indication du résultat final. mais bon.
ça ne joue pas, alors désolé de ne pas être d'accord avec les informations données. Difficiles de comprendre quand des personnes comme sacamerde te balance des truc en disant fait ça sans donner autre information.
Dommage mais bon tampis.

Merci snipeur pour ta compréhension et pour les explications, je ne dirait pas autant de sacamerde...

la moindre des chose et de répondre au questions et non pas de dire fais ça 4a marche et rien d'autre
0