le Rééchantillonnage des sons

 Pourquoi rééchantillonner ?

Les cartes sons PC utilisent des fréquences multiples de 11025 Hz (en particulier, 4*11025 = 44100 Hz du son dit « qualité CD »). Malheureusement, nos chères machines utilisent plutôt des fréquences multiples de 12500 Hz. Aussi, pour rejouer des sons issus du monde PC, il est nécessaire de procéder à ce que l’on appelle le rééchantillonnage des sons afin de faire comme si le son avait été enregistré à une fréquence multiple de 12500 Hz. Dans cet article, nous n’hésiterons pas à massacrer quelques aigus car nous allons voir comment procéder à ce rééchantillonnage par la méthode de l’interpollation linéaire. Et puisqu’on est sympa chez les Removers, en guise de dessert, nous vous dévoilerons la routine utilisée dans AviPlayer (celle dont je suis si fier dans la doc !)

 La méthode de l’interpollation linéaire

Imaginons donc avoir un son échantillonné à la fréquence (en Herz) f_r et devant être rejoué à la fréquence (toujours en Herz) f_p > f_r. La situation est la suivante : le son original est une belle courbe continue, le son échantillonné d’origine est la liste des valeurs de cette fonction apparaissant tous les 1/f_r secondes et le son rééchantillonné devrait correspondre à la liste des valeurs de cette fonction apparaissant tous les 1/f_p secondes. Le problème est : comment construire cette nouvelle liste de valeurs connaissant la première ?

L’une des méthodes les plus simples est celle que l’on appelle l’interpollation linéaire. L’idée est que pour reconstituer notre liste de valeurs, on fait comme si la courbe du son d’origine était la ligne brisée que définit le son échantillonné à la fréquence f_r. Pour obtenir un son échantillonné à la fréquence f_p, il n’y a donc plus qu’à lire les valeurs tous les 1/f_p secondes et le tour est joué (je vous avais bien dit que c’était simple).

 Calcul des nouveaux échantillons

Commençons par le calcul du n^e échantillon. Ce n^e échantillon « tombe » entre deux échantillons consécutifs du son original ; disons entre le k^e et le (k+1)^e échantillon. Soit v_k la valeur du k^e échantillon d’origine et v_k+1 celle de l’échantillon suivant. On se propose de calculer w_n, la valeur du n^e échantillon destination, en faisant comme s’il était sur le segment de droite défini par les deux échantillons d’origine. Si on veut formaliser un tout petit peu plus, on a donc un segment [A ;B] avec A=(k/fr ; v_k) et B=((k+1)/f_r ; v_k+1) et sur ce segment, on a un point C=(n/f_p ; w_n) dont on cherche à calculer l’ordonnée.

Pour réaliser ce calcul, on peut exploiter la particularité des vecteurs AC et AB : ils sont colinéaires. On a donc (n/f_p-k/f_r ; w_n - v_k) = t * ((k+1)/f_r - k/f_r ; v_k+1 - v_k), c’est-à-dire, après simplifications ((n*f_r -k*f_p) / (f_r * f_p) ; w_n - v_k) = t * (1/f_r ; v_k+1 - v_k). On en déduit que, après simplifications, t = n * (f_r / f_p) - k et donc finalement w_n = t*v_k+1 + (1-t)*v_k. A ce stade je sais que j’ai perdu 80% des lecteurs, mais rassurez-vous, tout va bien se passer.

Intéressons nous maintenant au calcul de k en fonction de n. Par définition, k est l’entier tel que k/f_r ≤ n/f_p < (k+1)/f_r, c’est-à-dire k est l’entier tel que k ≤ n * (f_r / f_p) < k+1. Autrement dit k est simplement la partie entière de n * (f_r / f_p). Si on injecte cette valeur dans t, on obtient que t = frac(n * (f_r / f_p)) (la partie fractionnaire de n * (f_r / f_p)).

 Implémentation pour des sons 8 bits signés

« Du code ! du code ! », je vous entends d’ici piaffer d’impatience... Bon, c’est bien parce que j’ai promis alors ! La première chose à remarquer est que le rapport f_r/f_p n’est pas entier (car f_p > f_r > 0) : il est même strictement compris entre 0 et 1. On va donc travailler en virgule fixe (on pourrait aussi travailler au copro en virgule flottante si on le voulait vraiment...) La virgule fixe, pour ceux qui ne le savent, consiste simplement à travailler avec des entiers multipliés par une puissance de 2 que l’on fixe (c’est très facile à manipuler ainsi avec des décalages). C’est à dire, le nombre décimal 0.5 est codé, avec une virgule fixe à 8 bits (c’est-à-dire en multipliant par 2^8 = 256), par 0.5*256 = 128. Evidemment, puisque l’on se donne une précision fixe (qui est la taille de la virgule), on perd de la précision. Par exemple 1/3 est codé par le quotient de 256 par 3, c’est-à-dire 85. La précision à laquelle on travaille est tout simplement donnée par 1/2^{(nombre de bits pour la virgule fixe)}.

Dans notre code, nous allons donc travailler avec une virgule fixe sur 8 bits. Le rapport f_r/f_p sera donc codé par le quotient de 256 * f_r / f_p.

Nous allons alors partir d’un registre de valeur nulle (nous sommes au début du son) et tant que nous n’aurons pas atteint la fin du son, nous allons ajouter d0 à ce registre. A chaque fois que ce registre dépassera la valeur d’un nouvel entier, nous avancerons d’un échantillon dans notre son original (cela correspond au fait que k est la partie entière de n * f_r/f_p) Si vous avez bien lu l’article de Stabylo sur les codes conditions, vous savez comment faire ceci en pratique. Il suffit de prendre un registre sur 8 bits, initialement à 0, et ajouter, sur 8 bits, d0 à ce registre. Si une retenue est générée (bit Carry positionné), cela signifie que nous avons dépassé un entier. Voilà donc le squelette de notre routine de rééchantillonnage :

Voilà, il ne nous reste plus qu’à compléter le code pour calculer le nouvel échantillon. Pour cela, nous avons en particulier besoin de la partie fractionnaire de n * f_r/f_p. Ca tombe bien, c’est justement ce que contiennent les 8 bits de poids faible de d1 : la partie fractionnaire est tout simplement ce qui est après la virgule ! Le code à compléter est donc :

Comment, quoi ? J’en entends au fond de la classe qui pouffent de rire et qui commencent à se moquer de moi. Comment ça « il y a plein de muls » ? Comment ça « ce code est lent » ? Comment ça « il n’y a pas de quoi être fier » ? Damned, je ne vous ai pas donné la version utilisé par AviPlayer qui a une boucle constituée — dixit la doc — de 11 instructions au lieu des 16 que je vous propose ici. Ok, ok, je vais voir ce que je peux faire alors... Disons que c’était pour commencer gentiment car après en fait, c’est que de la bidouille...

 Les précalculs : l’art d’être bourrin

Malheureusement, pas de secret, le boulot que l’on doit faire est bien celui décrit ci-dessus : le truc est qu’on peut en pré-mâcher un morceau avant de procéder effectivement au rééchantillonnage... Et ce n’est pas (trop) gourmand en mémoire car on a choisi une virgule fixe sur 8 bits (on aurait pu prendre un peu moins peut-être) et des échantillons sur 8 bits... En effet, comme t varie entre 0 et 1, il en est de même pour 1-t. Comme on travaille en virgule fixe sur 8 bits, cela veut dire que t et 1-t varient entre 0 et 255 (pour rappel, 2⁸=256). Il en est de même pour v_k (puisque les échantillons sont sur 8 bits).

L’idée de bourrin est donc, pour chaque valeur possible de t et 1-t et pour chaque valeur d’échantillon v possible, de pré-calculer t*v et (1-t)*v et de stocker ces résultats dans un tableau. Il suffit alors de modifier le code ci-dessus pour utiliser cette table plutôt que de faire ces calculs sur le moment. La dernière astuce est la lecture des valeurs dans la table de précalculs. La table stockée en mémoire représente donc un tableau à deux dimensions de la forme [valeurs pour t=0] [valeurs pour t=1] ... [valeurs pour t=256].

En fait, pour pré-mâcher un peu plus, on va intercaler dans cette table les valeurs calculées pour 1-t, c’est-à-dire [valeurs pour 1-t quand t=0] [valeurs pour 1-t quand t=1] ... [valeurs pour 1-t quand t=256] (on pourrait certes s’en passer puisque ces informations sont redondantes, mais ça permet de gagner encore un peu). La table est donc — pour fixer les idées — de la forme :

Chaque sous-tableau pèse en mémoire 2*256*16 bits (16 bits pour chaque valeur stockée), c’est à dire 1024=2¹⁰ octets ; ça nous fait donc un précalcul de 256 ko en tout (d’où le qualificatif de bourrin). Comme on ne veut pas trop perdre de temps à calculer l’endroit dans la table et qu’on a un 68020 (ou plus) à disposition, on va utiliser l’adressage fait exprès pour lire dans les tableaux. Le coefficient multiplicateur peut être de 2, 4 ou 8. Manque de chance, on aimerait multiplier par 1024. Heureusement, on peut aussi décider de stocker nos 8 bits additionneur non pas dans l’octet de poids faible du mot de poids faible des registres d0/d1 mais dans l’octet de poids fort de ce mot de poids faible. On gardera un comportement identique si l’on travaille alors sur des mots (en particulier pour la retenue et le bcc qui l’utilise). En faisant ça, on aura déjà prémultiplié par 2⁸ = 256 notre index. Et miracle, 1024 = 256*4 ! On peut donc effectivement utiliser le mode d’adressage 68020 et notre routine finale est la suivante :

On obtient donc une boucle principale de 11 instructions comme annoncé dans la doc de AviPlayer. A noter que quitte à faire des précalculs, on peut aussi traiter les problèmes de signe en construisant la table des valeurs (ce que j’ai supposé ici puisque j’ai retiré les ext.w)

 Conclusion

Nous avons vu dans cet article la méthode de rééchantillonnage par interpollation linéaire et vu comment l’implémenter en assembleur. Enfin, la routine d’AviPlayer a été disséquée et nous avons vu comment utiliser des précalculs.

 Addendum

Vous allez bien rigoler... hier, en jouant un peu avec le DSP de la Jaguar, j’ai trouvé une optimisation toute bête pour l’interpollation linéaire. En effet, il suffit de constater que w_n = t*v_k+1 + (1-t)*v_k = v_k + t * (v_k+1 - v_k). Et hop, une multiplication de gagnée ! Bien sûr, pour être valable, il faut travailler avec un bit de plus pour ne pas perdre des bits au passage. Dans le cas de la version avec précalculs, cela signifie que la moitié des précalculs ne sert plus à rien et ça économise un accès bus... pas mal non ?