El filtro de estado variable que describimos en la anterior entrega lo habíamos construido a partir de la discretización de la ecuación diferencial de un circuito electrónico. Ahora hablaremos de otro método para el diseño de filtros digitales que es bastante más intuitivo y mediante el cual podremos realizar implementaciones directas en digital.
La transformada z o H(z) de un sistema digital, es la ecuación que define totalmente su comportamiento. Si, por ejemplo, tenemos un sistema digital que responde a esta ecuación en diferencias: y[n] = 2·x[n] + x[n-1] - y[n-2], siendo x[n] la entrada en el instante n e y[n] la salida en el instante n, haciendo la transformada z de esta ecuación obtenemos: Y(z) = 2·X(z) + z-1·X(z) - z-2·Y(z). Como vemos: w[n-k] → z-k·W(z), siendo W(z) la transformada z de w[n] (la flecha → significa "su transformada z es") para cualquier w[n]. Vemos también que se trata de una transformación lineal, es decir: si tenemos I(z), J(z) y K(z), transformadas z de i[n], j[n] y k[n] respectivamente y tenemos que i[n] = a·j[n] + b·k[n], entonces I(z) = a·J(z) + b·K(z) (suponiendo a y b constantes).
Si en la ecuación Y(z) = 2·X(z) + z-1·X(z) - z-2·Y(z) dejamos a un lado de la igualdad la expresión Y(z)/X(z), en el otro lado de la igualdad obtendremos: (2 + z-1) / (1 + z-2). La ecuación de transferencia de un sistema discreto viene determinada por H(z) y ésta debe cumplir que Y(z)=X(z)·H(z) (por ser función de transferencia). Vemos que H(z)=Y(z)/X(z), por lo tanto la función de transferencia del sistema discreto es: H(z) = Y(z)/X(z) = (2 + z-1) / (1 + z-2).
En la mayoría de los casos la función de transferencia de un sistema discreto vendrá dada por una fracción de polinomios de diferente grado. En la función anterior podemos multiplicar el numerador y el denominador por z2 para mantener la igualdad, con lo que obtenemos H(z) = (2·z2 + z) / (z2 + 1) eliminando así los exponentes negativos.
Llamaremos "ceros" a los valores de z que hacen que H(z) = 0 y "polos" a los valores de z que hacen que H(z) = ∞. Los ceros serán las raices del polinomio del numerador mientras que los polos serán las raices del polinomio del denominador. En este caso los ceros son z=-1/2 y z=0 y los polos z=i y z=-i (z es una variable compleja).
En el plano complejo cada punto posee una parte real y una parte imaginaria (z = a + bi) y cualquier punto puede ser representado, además en forma de módulo y argumento, es decir, el punto z = a + bi, posee un módulo |z| = sqrt(a2 + b2) y un ángulo respectro del centro (z = 0) y del semieje horizontal positivo (argumento) arg(z) = arctan(b / a).
Cada uno de los ceros y polos del sistema se pueden localizar en el plano complejo y para cada uno de ellos se puede calcular su módulo y argumento.
Además de la ubicación de los ceros y los polos tenemos un elemento esencial dentro del plano complejo: la circunferencia goniométrica (una circunferencia con centro en el origen y de radio = 1). El comportamiento del sistema y el filtrado de frecuencias vendrá dado por la posición de los ceros y los polos con respecto a esta circunferencia.
Hay ciertas cosas a tener en cuenta: Los ceros o polos que posean una parte imaginaria diferente de 0 se encuentran por parejas de complejo conjugado (es decir, si tenemos un polo por ejemplo, en z = 1 + 2i también tendremos un polo en z = 1 - 2i, y lo mismo con los ceros) y esto es así debido a que las señales con las que trabajamos son señales sin componente imaginaria (señales reales), de esta forma los coeficientes de la ecuación en diferencias nos salen reales, sin parte imaginaria. Otra cosa a tener en cuenta es que cada cero o polo del sistema tiene asociada una frecuencia de señal en función de su argumento (el argumento de un cero o un polo estará comprendido entre 0 y PI, 0 representa la frecuencia 0 o componente de continua y PI representa la máxima frecuencia discreta o, lo que es lo mismo, la mitad de la frecuencia de muestreo).
Para el caso de los ceros, a medida que un cero se "aproxima" al borde de la circunferencia goniométrica desde el centro atenúa en mayor medida la frecuencia asociada. Un cero de módulo 1 (sobre la circunferencia) anula totalmente la frecuencia de ese cero (su argumento).
Para el caso de los polos, a medida que el polo se "aproxima" al borde de la circunferencia goniométrica desde el centro se produce una amplificación de la frecuencia asociada. En el caso de que un polo alcance la circunferencia (módulo = 1) se producirá una autooscilación del sistema a la menor excitación en la entrada (el sistema comenzará a oscilar a la frecuencia asociada al polo. ¡Podemos utilizar un filtro como oscilador digital!). Un polo que posea un módulo > 1 provocará un comportamiento inestable del sistema y oscilaciones cada vez más amplificadas en cada período.
Los ceros los podremos alojar en cualquier punto del plano complejo sin que ello provoque autooscilaciones o inestabilidades.
Recordemos los ceros y los polos del sistema anterior: Los ceros eran z=-1/2 y z=0 y los polos z=i y z=-i. Vemos que, al poseer un cero en z=-1/2 se producirá cierta atenuación de la componente de continua de la señal (los ceros en el origen son poco relevantes para nuestro caso), tenemos por otro lado dos polos complejos conjugados en i y -i (de módulo = 1), por lo tanto, el sistema autooscilará a la menor excitación en la entrada con una frecuencia PI/2, es decir, a una frecuencia f_max/2, siendo f_max la mitad de la frecuencia de muestreo.
El proceso es reversible, es decir, podemos partir de un diagrama de ceros y polos que diseñemos nosotros mismos y obtener la ecuación en diferencias del sistema realizando la transformada z inversa de la función H(z) = Y(z)/X(z). Podemos, pues, diseñar nuestros propios filtros y pasarlos automáticamente a ecuación en diferencias para ser implementados en el ordenador.
Otra cosa que hay que tener en cuenta a la hora de implementar un filtro digital es el escalado de H(z). Normalmente, cuando creamos un filtro dibujando los polos y los ceros sobre el plano z, luego hay que realizar un escalado de la función H(z) resultante para que, en la banda de paso (banda de frecuencias que dejamos pasar), H(z) = 1. El escalado de H(z) se realiza con una constante: Hfinal(z) = k · H(z), siendo k la constante calculada para que Hfinal(z) = 1 en la banda de paso del filtro.
Al lector interesado en estas cuestiones se le recomienda que acuda a bibliografía especializada (autores como Oggata y Manolakis) para profundizar. El espacio de esta serie de artículos no es suficiente para introducir todos los coneptos relacionados con el diseño de filtros digitales.
En el listado 1 podemos ver una implementación de este sistema discreto. Al ejecutarlo se puede comprobar como oscila a la frecuencia PI/2.
#include <stdio.h>
#include <stdlib.h>
#define TAM 500 /* generamos 500 muestras en total */
/* entrada en los instantes de tiempo n y n-1 */
/* para excitar el filtro la primera muestra la ponemos a 1.0 */
float input_n = 1.0, input_n_1 = 0;
/* salida en los instantes de tiempo n, n-1 y n-2 */
float output_n = 0, output_n_1 = 0, output_n_2 = 0;
int main(void) {
int i;
for (i = 0; i < TAM; i++) {
/* aplicamos la ecuación en diferencias */
output_n = 2 * input_n + input_n_1 - output_n_2;
/* mostramos la salida de forma numérica */
/* podremos ver cómo la onda de salida tiene una perioricidad de 4
muestras, es decir, una frecuencia de PI/2, la mitad de la
frecuencia máxima (la frecuencia máxima (PI) tiene una perioricidad
de 2 muestras)
si emitimos esta señal a través de la tarjeta de sonido con una
frecuencia de muestreo de 44.1 KHz, oiremos un tono de 11050 Hz */
printf("%.5f\n", output_n);
/* realizamos los desplazamientos de las muestras */
output_n_2 = output_n_1;
output_n_1 = output_n;
input_n_1 = input_n;
input_n = 0; /* anulamos las siguientes muestras de entrada */
}
return 0;
}
|
Listado 1 - Implementación de la ecuación en diferencias y[n] = 2·x[n] + x[n-1] - y[n-2]. |
Existen otros métodos de síntesis y de tipos de modulación además de los vistos en el artículo anterior. Si uno se fija en las publicaciones especializadas se dará cuenta de que los fabricantes de sintetizadores cada vez que sacan un nuevo producto al mercado lo venden como la panacea y la revolución sonora. Basta con ahondar un poco en las especificaciones para darse cuenta de que detrás de esas siglas incomprensibles se esconden ligeras variantes de métodos de síntesis triviales.
A continuación repasaremos algunos métodos de síntesis realmente diferentes de los vistos en el artículo anterior y que constituyen un resumen aproximado de lo que se puede encontrar hoy en el mundo de la síntesis.
La ecuación de generación de un tono digital como vimos en el artículo anterior se podía implementar así: f(n) = cos(2 · PI · F_SEÑAL · (n / F_MUESTREO)), siendo F_SEÑAL la frecuencia en Hz de la señal a reproducir, F_MUESTREO la frecuencia de muestreo del dispositivo y n el índice de la muestra a generar. Para explicar la distorsión de fase reescribamos este código: inicializamos una variable p a 0 (la fase) y ponemos la función f(n) = cos(2 · PI · p) en un bucle en n para generar cada una de las muestras. En cada iteración hacemos, además, p += Δp, (vamos incrementando la fase) siendo Δp = F_SEÑAL / F_MUESTREO. Por otro lado vigilamos que la variable p no se haga mayor que 1, haciendo que "de la vuelta" (ver el listado 2). De esta forma obtenemos también una señal senoidal de frecuencia F_SEÑAL a la frecuencia de muestreo F_MUESTREO.
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define F_SENIAL 220 /* un LA en la 3ª octava del piano */
#define F_MUEST 44100 /* frecuencia de muestreo: 44.1 KHz */
#define TIEMPO 2 /* segundos de reproducción */
#define MUESTRAS_TOTALES (TIEMPO * F_MUEST)
/* incremento en la variable fase, dicha variable va de 0 a 1 en cada período */
#define INC_FASE ((float)F_SENIAL / F_MUEST)
signed short int salida[MUESTRAS_TOTALES];
/* una distorsión sencillita */
float distorsion1(float fase) {
return (float)pow(fase, 3.0);
}
/* otra distorsión, esta vez definida como una función a trozos */
float distorsion2(float fase) {
if (fase < 0.5)
return (float)pow(fase, 2.0);
else if (fase < 0.75)
return (fase - 0.25);
else
return 0.80;
}
int main(void) {
int i;
float v;
float fase;
/* empezamos con fase 0 */
fase = 0;
for (i = 0; i < MUESTRAS_TOTALES; i++) {
/* calculamos la señal a partir de la fase aplicando un algoritmo de
distorsión que devuelva un valor entre 0 y 1 (como la variable fase,
que toma valores entre 0 y 1)
si no queremos que haya distorsión de fase simplemente sustituimos
la línea de código por v = cos(2 * M_PI * fase) */
v = cos(2 * M_PI * distorsion1(fase));
salida[i] = (signed short int) rint(v * 32767);
/* incrementamos la fase */
fase += INC_FASE;
/* mantenemos la variable fase dentro del intervalo [0, 1] */
if (fase > 1)
fase -= 1;
}
return 0;
}
|
Listado 2 - Ejemplo de síntesis mediante distorsión de fase. |
La distorsión de fase consiste en aplicar una función a la variable p anterior, de tal forma que f(n) = cos(2 · PI · DP(p)). DP(p) será una función cualquiera que, dependiendo del valor de su parámetro (p) devuelve un número entre 0 y 1 y que, en consecuencia, generará una distorsión en la fase de la señal de salida. Como vemos se parece mucho a una FM en la que hemos sustituido la señal moduladora por un algoritmo.
En el listado 2 tenemos un ejemplo de código que genera un tono mediante distorsión de fase. Se han incluido dos algoritmos de distorsión (distorsion1 y distorsion2).
La sincronización de onda consiste en dotar al oscilador de una entrada que, al ser activada, reinicia la fase de la onda a 0 (esta entrada se suele denominar sync). Una forma muy extendida es la de utilizar como señal sync la salida de otro oscilador a frecuencia audible y desplazada con respecto a la frecuencia del primer oscilador ([osc A] --sync--> [osc B]). De esta forma se generan componentes frecuenciales altas al introducir la sincronización saltos bruscos en la señal del oscilador sincronizado.
En el listado 3 tenemos un código que genera dos señales de onda cuadrada de frecuencias ligeramente diferentes (F_SENIAL y (F_SENIAL · MULT)) y que permite mezclarlas tal cual, o mezclarlas haciendo que la segunda de ellas se sincronice con la primera.
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define F_SENIAL 220 /* un LA en la 3ª octava del piano */
#define F_MUEST 44100 /* frecuencia de muestreo: 44.1 KHz */
#define MULT 1.008 /* factor de frecuencia para el segundo oscilador */
#define TIEMPO 2 /* segundos de reproducción */
#define MUESTRAS_TOTALES (TIEMPO * F_MUEST)
/* incrementos de la fase para cada oscilador */
#define INC_FASE1 ((float)F_SENIAL / F_MUEST)
#define INC_FASE2 ((float)(F_SENIAL * MULT) / F_MUEST)
#define SYNC 1 /* a 0 para no sincronizar, a 1 para sincronizar */
signed short int salida[MUESTRAS_TOTALES];
int main(void) {
int i;
float v1, v2, v;
float fase1, fase2;
/* inicializamos las fases a 0 */
fase1 = fase2 = 0;
for (i = 0; i < MUESTRAS_TOTALES; i++) {
/* primer oscilador */
v1 = cos(2 * M_PI * fase1);
v1 = (v1 > 0) ? 1.0 : -1.0;
/* segundo oscilador */
v2 = cos(2 * M_PI * fase2);
v2 = (v2 > 0) ? 1.0 : -1.0;
/* sumamos */
v = (v1 + v2) / 2.0;
salida[i] = (signed short int) rint(v * 32767);
/* incrementamos la fases y sincronizados, si es preciso */
fase2 += INC_FASE2;
if (fase2 > 1)
fase2 -= 1;
fase1 += INC_FASE1;
if (fase1 > 1) {
fase1 -= 1;
if (SYNC)
fase2 = 0;
}
}
return 0;
}
|
Listado 3 - Ejemplo de síntesis mediante sincronización de onda. |
Modelar una onda consiste en utilizar una función g y aplicarla a la señal de salida. Si tenemos que la señal de salida de un oscilador, de un conjunto de osciladores, o, en general, de un sistema de síntesis es x[n], podemos hacer que la señal de salida y[n] = g(x[n]), siendo g(x) una función cualquiera. Si g(x) = x, tenemos que y[n] = x[n]. Imaginemos la función g(x) = 1 si x > 0, g(x) = -1 si x i< 0. Si realizamos el modelado de onda y[n] = g(x[n]) y x[n] es una onda senoidal, y[n] será la versión cuadrada de x[n]. Las posibilidades son infinitas, pudiendo, por ejemplo, utilizar como g(x) una arcotangente, una spline o lo que nos apetezca.
La síntesis vectorial es un método muy sencillo y a la vez muy versátil de síntesis basado en la mezcla de varias formas de onda de forma lineal. Para comprender este método de síntesis, supongamos que tenemos un cuadrado y un punto dentro de ese cuadrado p y a cada una de las esquinas del cuadrado le asignamos un oscilador con una forma de onda diferente. La forma de onda resultante será la suma ponderada de las cuatro formas de onda de las esquinas en función de la posición del punto p. Si, por ejemplo, colocamos p en la esquina correspondiente al oscilador A, sólo oiremos como salida la onda del oscilador A. Situando el punto p en cualquiera de los puntos centrales del cuadrado obtendremos diferentes timbres en función de la distancia entre p y cada uno de los osciladores (a menor distancia, mayor presencia de ese oscilador en la señal de salida).
Este método de síntesis se basa en la utilización de muchas formas de onda diferentes controlando la mezcla final mediante LFOs o envolventes. Imaginemos una escala del 0 al 1. En la posición 0 ponemos un oscilador de diente de sierra, en la posición 0.33 un oscilador de onda cuadrada y en la posición 1 ponemos un oscilador de onda triangular. Las posiciones intermedias, al igual que en la síntesis vectorial, dan como resultado una mezcla lineal de las tres formas de onda en función de la distancia entre el punto de la escala donde nos encontremos y cada oscilador.
Si hacemos que una envolvente ADRS o AD, o un LFO haga variar este punto entre 0 y 1 obtendremos un sonido que variará entre las diferentes formas de onda. Las formas de onda pueden ser del tipo que queramos.
En el listado 4 tenemos un ejemplo de código que genera un sonido mediante secuenciación de onda entre 3 formas de onda.
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define F_SENIAL 220 /* LA en la 3ª octava del piano */
#define F_MUEST 44100 /* frecuencia de muestreo */
#define TIEMPO 2 /* duración total en segundos */
#define TIEMPO_ATAQUE 0.5 /* tiempo de ataque en segundos */
#define TIEMPO_CAIDA 1.0 /* tiempo de caida en segundos */
#define MUESTRAS_TOTALES (TIEMPO * F_MUEST)
#define MUESTRAS_ATAQUE ((int) rint(TIEMPO_ATAQUE * F_MUEST))
#define MUESTRAS_CAIDA ((int) rint(TIEMPO_CAIDA * F_MUEST))
/* variables para el control del oscilador de diente de sierra */
int ds_t = 0;
int ds_max_t;
float ds_ramp;
/* control de la envolvente */
int env_t = 0;
/* arrays donde meteremos las diferentes señales */
float salida_diente_sierra[MUESTRAS_TOTALES];
float salida_cuadrada1[MUESTRAS_TOTALES];
float salida_cuadrada2[MUESTRAS_TOTALES];
float salida_envolvente[MUESTRAS_TOTALES];
float mezcla[MUESTRAS_TOTALES];
/* array de salida para enviar directamente al dispositivo */
signed short int salida[MUESTRAS_TOTALES];
/******************************************************************************
Inicializa las variables para generar la onda de diente de sierra en tiempo
real mediante el cálculo de la rampa.
*/
void diente_sierra_inicia(float frec) {
ds_t = 0;
/* ds_max_t es el número de muestras por período de la señal */
ds_max_t = (int) rint(F_MUEST / frec);
ds_max_t--;
/* pendiente del diente de sierra */
ds_ramp = 2.0 / ds_max_t;
return;
}
/******************************************************************************
Genera n muestras en el array out de la señal de diente de sierra.
*/
void diente_sierra_trabaja(float *out, int n) {
int i;
for (i = 0; i < n; i++) {
/* generamos la salida multiplicando la t del oscilador por la rampa */
out[i] = -1.0 + (ds_t * ds_ramp);
ds_t++;
if (ds_t > ds_max_t)
ds_t = 0;
}
return;
}
/******************************************************************************
Genera, a partir de una señal de entrada, una onda cuadrada con anchura
ajustable a pwm.
*/
void cuadrada_genera(float *out, float *ds_in, int n, float pwm) {
int i;
for (i = 0; i < n; i++)
/* si la señal de entrada supera pwm emitimos un +1, si no un -1 */
out[i] = (ds_in[i] > pwm) ? 1.0 : -1.0;
return;
}
/******************************************************************************
Genera n muestras sobre el array out de una envolvente de tipo AD (valores
de salida entre 0 y 1).
*/
void envolvente_genera(float *out, int n) {
int i;
for (i = 0; i < n; i++) {
/* período de ataque */
if (env_t <= MUESTRAS_ATAQUE)
out[i] = (float)env_t / MUESTRAS_ATAQUE;
/* período de caida */
else if (env_t <= (MUESTRAS_ATAQUE + MUESTRAS_CAIDA)) {
out[i] = ((float)(env_t - MUESTRAS_ATAQUE) / MUESTRAS_CAIDA);
out[i] = 1.0 - out[i];
}
/* fin de la envolvente */
else
out[i] = 0;
env_t++;
}
return;
}
/*****************************************************************************/
int main(void) {
int i;
float nivel1, nivel2, nivel3;
/* generamos las diferentes señales sobre los arrays */
diente_sierra_inicia(F_SENIAL);
diente_sierra_trabaja(salida_diente_sierra, MUESTRAS_TOTALES);
cuadrada_genera(salida_cuadrada1, salida_diente_sierra, MUESTRAS_TOTALES,
0);
cuadrada_genera(salida_cuadrada2, salida_diente_sierra, MUESTRAS_TOTALES,
0.8);
envolvente_genera(salida_envolvente, MUESTRAS_TOTALES);
/* hacemos la secuenciación en función de la envolvente */
for (i = 0; i < MUESTRAS_TOTALES; i++) {
/* expandimos la envolvente de [0, 1] a [-1, +1] */
salida_envolvente[i] = (2.0 * salida_envolvente[i]) - 1.0;
/* calculamos los niveles para cada onda */
if (salida_envolvente[i] < 0) {
nivel3 = 0.0;
nivel1 = -salida_envolvente[i];
nivel2 = 1.0 - nivel1;
}
else {
nivel1 = 0.0;
nivel3 = salida_envolvente[i];
nivel2 = 1.0 - nivel3;
}
/* calculamos la mezcla final de señales */
mezcla[i] = (nivel1 * salida_diente_sierra[i]) +
(nivel2 * salida_cuadrada1[i]) +
(nivel3 * salida_cuadrada2[i]);
/* convertimos la mezcla final para emitirla por el dispositivo */
if (mezcla[i] > 1.0)
mezcla[i] = 1.0;
else if (mezcla[i] < -1.0)
mezcla[i] = -1.0;
salida[i] = (int) rint(mezcla[i] * 32767);
}
return 0;
}
|
Listado 4 - Ejemplo de síntesis mediante secuenciación de onda (también llamado síntesis vectorial) utilizando una envolvente. |
En el crossfade a frecuencia audible lo que hacemos simplemente es una secuenciación de onda mediante un LFO que oscila a frecuencia audible y entre 2 formas de onda situadas en los puntos extremos (0 y 1) de la escala (en este caso hablar de LFO no sería correcto, ya que genera frecuencias audibles).
Tradicionalmente, en los sintetizadores analógicos, la forma de onda que gobierna el crossfade suele ser una onda cuadrada con diferentes anchuras de pulso. Utilizando una onda cuadrada lo que obtenemos a la salida es o la señal de un oscilador o la señal del otro, sin mezclas. Una forma más versátil, es utilizar formas de onda más extrañas, como una triangular, un coseno rectificado a onda completa o a media onda, etc.
En el listado 5 tenemos un programa que genera un crossfade a frecuencia audible entre dos formas de onda sencillas.
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define F_SENIAL 220 /* LA en la 3ª octava del piano */
#define F_CROSSFADE 111 /* frecuencia del crossfade */
#define F_MUEST 44100 /* frecuencia de muestreo */
#define TIEMPO 2 /* duración total en segundos */
#define TIEMPO_ATAQUE 0.5 /* tiempo de ataque en segundos */
#define TIEMPO_CAIDA 1.0 /* tiempo de caida en segundos */
#define MUESTRAS_TOTALES (TIEMPO * F_MUEST)
/* variables para el control del oscilador de diente de sierra */
int ds1_t = 0;
int ds1_max_t;
float ds1_ramp;
int ds2_t = 0;
int ds2_max_t;
float ds2_ramp;
/* arrays donde meteremos las diferentes señales */
float salida_diente_sierra[MUESTRAS_TOTALES];
float salida_cuadrada[MUESTRAS_TOTALES];
float salida_crossfade[MUESTRAS_TOTALES];
float mezcla[MUESTRAS_TOTALES];
/* array de salida para enviar directamente al dispositivo */
signed short int salida[MUESTRAS_TOTALES];
/******************************************************************************
Inicializa las variables para generar la onda de diente de sierra en tiempo
real mediante el cálculo de la rampa.
*/
void diente_sierra1_inicia(float frec) {
ds1_t = 0;
/* ds_max_t es el número de muestras por período de la señal */
ds1_max_t = (int) rint(F_MUEST / frec);
ds1_max_t--;
/* pendiente del diente de sierra */
ds1_ramp = 2.0 / ds1_max_t;
return;
}
void diente_sierra2_inicia(float frec) {
ds2_t = 0;
/* ds_max_t es el número de muestras por período de la señal */
ds2_max_t = (int) rint(F_MUEST / frec);
ds2_max_t--;
/* pendiente del diente de sierra */
ds2_ramp = 2.0 / ds2_max_t;
return;
}
/******************************************************************************
Genera n muestras en el array out de la señal de diente de sierra.
*/
void diente_sierra1_trabaja(float *out, int n) {
int i;
for (i = 0; i < n; i++) {
/* generamos la salida multiplicando la t del oscilador por la rampa */
out[i] = -1.0 + (ds1_t * ds1_ramp);
ds1_t++;
if (ds1_t > ds1_max_t)
ds1_t = 0;
}
return;
}
void diente_sierra2_trabaja(float *out, int n) {
int i;
for (i = 0; i < n; i++) {
/* generamos la salida multiplicando la t del oscilador por la rampa */
out[i] = -1.0 + (ds2_t * ds2_ramp);
ds2_t++;
if (ds2_t > ds2_max_t)
ds2_t = 0;
}
return;
}
/******************************************************************************
Genera, a partir de una señal de entrada, una onda cuadrada con anchura
ajustable a pwm.
*/
void cuadrada_genera(float *out, float *ds_in, int n, float pwm) {
int i;
for (i = 0; i < n; i++)
/* si la señal de entrada supera pwm emitimos un +1, si no un -1 */
out[i] = (ds_in[i] > pwm) ? 1.0 : -1.0;
return;
}
/*****************************************************************************/
int main(void) {
int i;
float nivel1, nivel2;
/* generamos las diferentes señales sobre los arrays */
diente_sierra1_inicia(F_SENIAL);
diente_sierra1_trabaja(salida_diente_sierra, MUESTRAS_TOTALES);
cuadrada_genera(salida_cuadrada, salida_diente_sierra, MUESTRAS_TOTALES,
0.7);
diente_sierra2_inicia(F_CROSSFADE);
diente_sierra2_trabaja(salida_crossfade, MUESTRAS_TOTALES);
/* hacemos el crossfade */
for (i = 0; i < MUESTRAS_TOTALES; i++) {
/* comprimimos el rango del crossfade de [-1, +1] a [0, 1] */
salida_crossfade[i] = (salida_crossfade[i] + 1.0) / 2.0;
/* calculamos los niveles para cada onda */
nivel1 = 1.0 - salida_crossfade[i];
nivel2 = salida_crossfade[i];
/* calculamos la mezcla final de señales */
mezcla[i] = (nivel1 * salida_diente_sierra[i]) +
(nivel2 * salida_cuadrada[i]);
/* convertimos la mezcla final para emitirla por el dispositivo */
if (mezcla[i] > 1.0)
mezcla[i] = 1.0;
else if (mezcla[i] < -1.0)
mezcla[i] = -1.0;
salida[i] = (int) rint(mezcla[i] * 32767);
}
return 0;
}
|
Listado 5 - Crossfade a frecuencia audible. |
La síntesis mediante ondas de terreno es un método de síntesis algo engorroso y con una complejidad de cálculo bastante grande comparado con el resto de métodos de síntesis que hemos visto.
Definimos matemáticamente una función f(x, y) que, para dos valores (x e y, semejante a lo que serían la latitud y la longitud en coordenadas terrestres) nos da otro valor (que podríamos llamar "altura") que se utiliza como salida de audio.
Los valores x e y son, a su vez señales de audio provenientes de un oscilador o de una tabla de ondas que van variando con el tiempo (x(t) e y(t)). Si suponemos que x(t)=A·cos(wt) y que y(t)=A·sen(wt), la trayectoria de los puntos (x, y) a lo largo del tiempo será la de una circunferencia de radio A y el punto (x, y) girará con una velocidad angular de w radianes/segundo.
Aplicando x(t) e y(t) a la función de 2 variables f(x, y) obtendremos un sonido que dependerá de la ecuación de f y de las funciones x(t) e y(t). Si aplicamos al valor de A una envolvente AD que conseguiremos que el punto (x, y), en lugar de describir una circunferencia todo el tiempo, parta del centro (A=0) y describa una espiral hasta A=1 (ataque) para luego volver a describir una espiral hasta el centro (caida).
En el listado 6 tenemos un código que genera un tono a partir de la ecuación g(x, y) = (x - y)·(x - 1)·(x + 1)·(y - 1)·(y + 1).
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define F_XY 220 /* un LA en la 3ª octava del piano */
#define F_MUEST 44100 /* frecuencia de muetreo: 44.1 KHz */
#define TIEMPO 2 /* segundos de reproducción */
#define ATAQUE 0.04
#define CAIDA 0.9
#define MUESTRAS_TOTALES (TIEMPO * F_MUEST)
#define MUESTRAS_ATAQUE ((int)(ATAQUE * F_MUEST))
#define MUESTRAS_CAIDA ((int)(CAIDA * F_MUEST))
signed short int salida[MUESTRAS_TOTALES];
int env_t = 0;
float envolvente(void) {
float v;
v = 0;
if (env_t <= MUESTRAS_ATAQUE) {
v = (float)env_t / MUESTRAS_ATAQUE;
}
else if ((env_t > MUESTRAS_ATAQUE) &&
((env_t - MUESTRAS_ATAQUE) <= MUESTRAS_CAIDA)) {
v = 1.0 - ((float)(env_t - MUESTRAS_ATAQUE) / MUESTRAS_CAIDA);
}
env_t++;
return v;
}
float g1(float x, float y) {
return (sqrt(x + 1) * y * x * (x - 4) * (y + 2));
}
float g2(float x, float y) {
return ((x - y) * (x - 1) * (x + 1) * (y - 1) * (y + 1));
}
int main(void) {
int i;
float v, x, y, env;
for (i = 0; i < MUESTRAS_TOTALES; i++) {
env = envolvente();
x = env * cos(2 * M_PI * ((float)i / F_MUEST) * F_XY);
y = env * sin(2 * M_PI * ((float)i / F_MUEST) * F_XY);
v = g2(x, y);
v = (v > 1.0) ? 1.0 : ((v < -1.0) ? -1.0 : v);
salida[i] = (int) rint(v * 32767);
}
return 0;
}
|
Listado 6 - Ejemplo de síntesis mediante ondas de terreno. |
La síntesis mediante tablas de ondas es el método de síntesis que se utiliza en los samplers y en los programas trackers muy utilizamos en ordenadores. Se basa en la reproducción de sonidos previamente grabados; ya sean de instrumentos reales, de ruidos o de cualquier tipo de sonido.
El principio es muy simple. Si, por ejemplo, grabamos una nota Do4 de una flauta travesera a una frecuenca de muestreo de 44100 Hz. Podremos obtener, al menos en teoría, toda la escala de notas del instrumento reproduciendo luego esas muestras a diferentes velocidades (a 22050 Hz si queremos oir un Do3 o a 88200 Hz si queremos oir un Do5, y así con el resto de notas de la escala). Para obtener la frecuencia en muestras/segundo a la que debemos recorrer un sonido muestreado: Veloc = F_MUEST_ORIG · pow(2, ((Nota - NOTA_ORIG) / 12) + (Octava - OCTAVA_ORIG)). Siendo Veloc la nueva frecuencia de reproducción que queremos calcular en muestras/segundo, F_MUEST_ORIG la frecuencia a la que fue muestreado el sonido, NOTA_ORIG la nota musical original del instrumento digitalizado (0=Do ... 11=Si), Nota la nota musical que queremos reproducir, OCTAVA_ORIG la octava en la que fue reproducida la nota original y Octava la octava a la que queremos reproducir la nota.
Como se puede apreciar, la ecuación anterior nos puede dar un valor por encima de la frecuencia de muestreo del dispositivo. Para evitarnos este tipo de problemas, es mejor trabajar utilizando incrementos en lugar de velocidades. Es más útil calcular el salto entre muestras que hay que realizar (ver el artículo anterior de la serie).
Para calcular el salto o incremento entre muestras, basta con dividir el resultado de la ecuación anterior (Veloc) entre la frecuencia de muestreo del dispositivo de salida (F_MUEST): Inc = Veloc / F_MUEST. Inc será, por lo general, un número fraccionario.
De esta forma, para reproducir las muestras de un instrumento previamente muestreado, calcularemos el valor de su incremento (Inc) para la nota que queramos reproducir, inicializamos una variable t = 0 que nos irá indicando el índice de la tabla y, en cada iteración del bucle de reproducción, emitimos el valor tabla[t] y hacemos t += Inc. Evidentemente, es recomendable utilizar algún método de interpolación para calcular el valor de tabla[t]. En el array tabla tendremos almacenadas las muestras del instrumento.
Normalmente, en la síntesis mediante tablas de ondas se introducen bucles en la onda, que permiten crear sonidos sostenidos y también se aplican envolventes de amplitud, de panorámica y de filtro a la onda. También se usan LFOs para vibrato (amplicándolos al valor de Inc) y para trémolo (aplicándolos a la amplitud).
En teoría, partiendo de una sola onda de un instrumento real, y variando la velocidad de recorrido, obtendríamos toda la gama de notas que es capaz de reproducir ese instrumento. Esto no es del todo cierto debido al proceso natural del sonido.
Cuando analizamos el sonido de un instrumento musical, este está compuesto por sus parciales (fundamental, 2·fundamental, 3·fundamental, 8·fundamental, etc) tanto armónicos (factor entero) como inarmónicos (factor no entero) y por los formantes, que son los parciales predominantes y que dependen de las resonancias propias del instrumento. Los formantes de un instrumento acústico se mantienen fijos (dependen de las características físicas del instrumento) y, a medida que reproducimos unas notas y otras, los parciales que coincidan con los formantes se amplificarán y los otros se atenuarán. En el caso de que hagamos simplemente una variación de la velocidad de recorrido de una onda muestreada, desplazaremos por igual parciales y formantes, provocando, por ejemplo, que un Do4 muestreado de un piano y posteriormente reproducido al doble de velocidad poco o nada tenga que ver con el sonido de un Do5 del mismo piano.
El multimuestreo surge como medida paliativa de este problema. El truco está, simplemente en muestrear varias notas de un mismo instrumento distribuyéndolas de forma equidistante a lo largo de toda la escala. En el caso de un piano podríamos muestrear 6 notas: Fa1, Fa2, Fa3, Fa4, Fa5 y Fa6, y obtener el resto de notas del piano en función del segmento donde estén alojadas. Por ejemplo, un Re5 lo obtendremos desplazando 3 semitonos hacia abajo el Fa5, mientras que el Si4 lo obtendríamos desplazando 6 semitonos hacia arriba (al agudo) el Fa4.
Entre los códigos fuente de este artículo, el fichero readxi.c es un programa que lee ficheros en formato .XI (eXtended Instrument, formato muy extendido del software FastTracker basado en tablas de ondas mediante multimuestreo y envolventes) y los reproduce a la nota que queramos.
La síntesis Karplus-Strong es uno de los primeros antecedentes de la síntesis de modelado físico. Se utiliza principalmente para la síntesis de sonidos percusivos o de cuerda pulsada y está basada en el uso de una línea de retardo realimentada. Si hacemos que la línea de retardo tenga una longitud de T segundos (T·F_MUEST muestras), conseguiremos un sonido con una frecuencia fundamental de 1/T Hz, las muestras que salen por un lado de la línea de retardo pasan por un filtro (en general, un H(z)) antes de entrar por el otro extremo. Si modelamos el filtro H(z) como un filtro paso-bajo con una ganancia en la banda de paso un poco menor que 1 e inicializamos la línea de retardo con una onda cuadrada o de pulso, lo que obtendremos será un sonido pulsado o percutido que evoluciona de una onda cuadrada y que, debido al filtro en la realimentación, va perdiendo las frecuencias altas y convirtiéndose, poco a poco en una onda senoidal que, finalmente, se extingue (al ser la ganancia en la banda de paso menor que 1).
Otra posibilidad es la de inicializar la línea de retardo con ruido, de esta forma oiremos como el ruido va perdiendo componentes de alta frecuencia hasta convertirse en un tono armónico. Como vemos las posibilidades son casi infinitas, ya que en H(z) podemos modelar cualquier tipo de filtro (paso-bajo, paso-alto, paso-banda, elimina-banda o cualquier otro) para obtener el espectro que queramos y, a su vez, la línea de retardo la podemos inicializar con ondas de las formas que queramos (cuadradas, triangular, ruido, etc).
Nos pararemos ahora un poco en la síntesis de sonidos de percusión. Existen, básicamente, 3 formas diferentes de sintetizar sonidos de percusión. La forma más sencilla sonsiste en muestrearlos y reproducirlos luego a la misma frecuencia de muestreo (este método no se puede considerar de síntesis). Los otros métodos son el modelado físico y la síntesis analógica. Nos centraremos en la síntesis analógica, al ser este tipo de síntesis el más utilizado junto con el muestreo.
Los bombos sintetizados tienen todos un denominador común: consisten en una onda senoidal cuyo tono cae rápidamente, de frecuencias del orden de los cientos de hercios hasta frecuencias del orden de los 30 o 40 hercios en un breve espacio de tiempo (milisegundos). Normalmente la caida es exponencial (a medida que va bajando la frecuencia también va bajando la velocidad de caida).
Elementos que se suelen añadir son: Un tick inicial de alta frecuencia para emular el sonido de un bombo real o la utilización, en el momento inicial de una onda no senoidal, con bastante contenido armónico. A veces también se utiliza un modelador de onda para distorsionar la forma senoidal y añadirle ese contenido armónico adicional.
Las cajas analógicas se sintetizan mediante la mezcla de una o más ondas senoidales (normalmente desafinadas entre sí y con una ligera caida tonal; más ligera que en el bombo), junto con ruido de alta frecuencia (ruido blanco filtrado mediante un paso-alto o un paso-banda, normalmente).
Los sonidos metálicos se obtienen mezclando ruido de alta frecuencia con otras ondas que aporten textura metálica al sonidos. Para aportar este tipo de textura se suele recurrir a la modulación en anillo y al filtrado de ondas pulso o cuadradas mezcladas entre sí. Otro instrumento muy curioso son las palmadas, que se sintetizan mediante el disparo en un corto espacio de tiempo (milisegundos) de 3 o más ruidos percusivos de alta frecuencia.
Entre los códigos de ejemplo hay un programa que se llama dsynth.c y que implementa estos tópicos para generar un ritmo estilo tecno.
La síntesis mediante modelado físico es un método basado en modelos de comportamiento de instrumentos reales ya sea para emularlos o para crear nuevos timbres.
Un método directo consiste en analizar el comportamiento de un instrumento real y describirlo mediante ecuaciones diferenciales. Luego basta con discretizar esas ecuaciones diferenciales y ya podremos emular el sonido de ese instrumento. Desgraciadamente la cosa no es tan fácil como parece. Normalmente las ecuaciones diferenciales hasta del instrumento más sencillo son un autentico calvario tanto para entenderlas como para implementarlas al poseer gran cantidad de parámetros (longitud y grosor de las cuerda, grosor de la madera, coeficientes relacionados con las características físicas de cajas de resonancia, interacción entre cuerdas, etc).
Una buena aproximación a la síntesis mediante modelado físico directo es la síntesis mediante guías de ondas (o waveguides, como se las conoce). Una guía de ondas es un registro de desplazamiento de muestras que simboliza la trayectoria de una onda. Veámoslo con un ejemplo.
Pensemos en una cuerda de guitarra. Se trata de una cuerda con dos puntos de apoyo en sus extremos y en la que, como sabemos, las ondas se desplazan en ambos sentidos. Cada sentido de propagación lo representaremos mediante un registro de desplazamiento (un array al fin y al cabo) cuya longitud lógica (en muestras) se corresponderá con la longitud física de la cuerda que queremos modelar (tendremos, pues, dos registros de desplazamiento que irán desplazando muestras en sentidos opuestos). Por otro lado tenemos los apoyos de la cuerda situados en los extremos de ésta. En los puntos de apoyo las muestras que salen de un registro de desplazamiento se meten por el otro registro de desplazamiento (recordemos que desplazan las muestras en sentidos opuestos), esto sucede en ambos extremos de la cuerda y, además, en los extremos, la cuerda no vibra.
Sean reg1[] y reg2[] los registros de desplazamiento que representan la cuerda (de tamaño n cada uno). En cada iteración del bucle de generación de sonido desplazamos, en reg1, las muestras hacia la derecha una posición y, en reg2, hacia la izquierda una posición. El valor que sale por la derecha del registro de desplazamiento reg1 lo introducimos, cambiado de signo, al final de reg2 (en reg2[n-1]) y el valor que sale por la izquierda del registro de desplazamiento reg2 lo instroducimos, cambiado de signo, en reg1[0]. La salida de audio la tomaremos de un punto k con 0≤k<n cualquiera siendo Salida = reg1[k] + reg2[k]. ¿Qué conseguimos con esto? Vemos que, mediante el cambio de signo que hacemos cuando una muestra sale de un registro para meterse en el otro, mantenemos las condiciones de vibración de una cuerda real, ya que, si k=0 ó k=n-1 (en los extremos) Salida=0 (es decir, no se emite sonido en los puntos extremos de la cuerda).
Para arrancar el modelo basta con inicializar los registros de desplazamiento con valores que correspondan a las condiciones iniciales de una cuerda real. Si, por ejemplo, queremos emular una cuerda pulsada con una pua inicializaremos las muestras dentro de los registros de desplazamiento en forma de rampa situando el pico en punto de la cuerda en el que actuamos con la pua.
Evidentemente, este modelo que acabamos de describir corresponde a el de una cuerda ideal, en una cuerda real, se produce una caida paulatina tanto de la amplitud como del contenido armónico del sonido inicial. Para provocar pérdidas en la amplitud lo que se hace es hacer que las muestras, al pasar de un registro de desplazamiento a otro (en los extremos), se atenúen un poco en cada vuelta y para modelar el contenido armónico lo que se hace es introducir filtros (normalmente, paso bajo) en los extremos de la cuerda.
En el listado 7 tenemos un código que genera un sonido de cuerda pulsada. En uno de los extremos de la cuerda hemos colocado un filtro H(z) con un cero en z=-1 como filtro sencillo paso bajo.
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
/* tamaño del buffer */
#define BUFFER_SIZE 256
/* frecuencia de muestreo */
#define SAMPLING_RATE 44100
/* frecuencia audible más baja que podremos reproducir */
#define LOWEST_MUSIC_FREQ 20.0
/* calculamos el tamaño máximo de la guía de ondas a partir de la frecuencia
audible más baja (LOWEST_MUSIC_FREQ) */
#define MAX_WAVEGUIDE_SIZE ((int)((1.0/LOWEST_MUSIC_FREQ) * SAMPLING_RATE/2.0))
typedef struct {
float freq;
int waveguide_size;
int mic_position;
int pluck_position;
float pluck_depth;
float forward[MAX_WAVEGUIDE_SIZE];
float backward[MAX_WAVEGUIDE_SIZE];
int start;
float filter_out;
float attenuation;
} plucked_string;
/******************************************************************************
Inicializa los campos de la estructura plucked_string para comenzar la
reproducción.
Entrada:
freq = frecuencia musical
mic_pos = posición normalizada (0..1) de la pastilla/micrófono
pluck_position = posición normalizada (0..1) en la que actúa la púa
pluck_depth = amplitud de la deformación inicial de la púa
attenuation = nivel de atenuación causado por el puente
*/
void plucked_string_init(plucked_string *p, float freq, float mic_pos,
float pluck_pos, float pluck_depth,float attenuation) {
/* calculamos el tamaño que deberá tener la guía de ondas para esta
frecuencia */
p->waveguide_size = (int) rint((1.0 / freq) * SAMPLING_RATE / 2.0);
/* denormalizamos las posiciones de la pastilla/micrófono y de la púa */
p->mic_position = (int) rint(mic_pos * (p->waveguide_size - 1));
p->pluck_position = (int) rint(pluck_pos * (p->waveguide_size - 1));
/* otras variables */
p->pluck_depth = pluck_depth;
p->filter_out = 0;
p->start = 1;
p->attenuation = attenuation;
return;
}
/******************************************************************************
Calcula BUFFER_SIZE muestras de la cuerda pulsada y las escribe en el
buffer output.
*/
void plucked_string_work(plucked_string *p, float *output) {
int i, k;
float forward_out, backward_out;
for (k = 0; k < BUFFER_SIZE; k++) {
if (p->start) {
/* dibujamos un triángulo en el array forward para simular la
pulsación de la cuerda */
/* el vértice superior del triángulo estará en p->pluck_position y
tendrá una altura de p->pluck_depth */
for (i = 0; i < p->pluck_position; i++)
p->forward[i] = ( ((float) i) / p->pluck_position ) *
p->pluck_depth;
for (i = p->pluck_position; i < p->waveguide_size; i++)
p->forward[i] = ( (((float) i) - p->waveguide_size) /
(p->pluck_position - p->waveguide_size) )
* p->pluck_depth;
/* en el array backward dibujamos la señal inversa del forward */
/* se hace esto para que en el instante inicial salida=0 y evitar
un salto brusco de la señal en el ataque */
for (i = 0; i < p->waveguide_size; i++)
p->backward[i] = -p->forward[i];
p->filter_out = 0;
p->start = 0;
}
else {
/* desplazamos los arrays */
forward_out = p->forward[p->waveguide_size - 1];
memmove(&(p->forward[1]), &(p->forward[0]),
(p->waveguide_size - 1) * sizeof(float));
backward_out = p->backward[0];
memmove(&(p->backward[0]), &(p->backward[1]),
(p->waveguide_size - 1) * sizeof(float));
/* en un extremo suponemos un apoyo ideal y reflejamos toda la
energía */
p->backward[p->waveguide_size - 1] = -forward_out;
/* ...en el otro extremo implementamos un filtro digital con un cero
en z=0 y un polo en z=1/2 (paso-bajo) y aplicamos una pequeña
atenuación (p->attenuation) */
p->filter_out = (p->filter_out - backward_out) / 2;
p->forward[0] = p->filter_out * p->attenuation;
}
/* la salida es la suma de las oscilaciones hacia uno y otro lado en un
mismo punto (p->mic_position) */
output[k] = p->forward[p->mic_position] + p->backward[p->mic_position];
output[k] /= 2;
}
return;
}
/*****************************************************************************/
int main(void) {
int i;
plucked_string ps;
float s[BUFFER_SIZE];
/* preparamos la cuerda pulsada */
plucked_string_init(&ps, 90, 0.2, 0.4, 1, 0.99);
/* generamos las muestras */
for (i = 0; i < 2000; i++) {
plucked_string_work(&ps, s);
/* en "s" tenemos BUFFER_SIZE muestras que podemos enviar al dispositivo de salida */
}
return 0;
}
|
Listado 7 - Ejemplo de síntesis por modelado físico mediante guías de ondas. Simulación de una cuerda pulsada con una púa. |
En esta segunda entrega hemos afianzado los conceptos sobre filtros digitales y nos hemos dado una vuelta por los métodos de síntesis menos conocidos, pero que aportan texturas interesantes al sonido. También hemos visto la percusión analógica; muy utilizada en la música actual y hemos hablado del modelado físico mediante guías de ondas y de las posiblidades que permite en el campo de la emulación de instrumentos reales.
En la próxima y última entrega hablaremos de la síntesis mediante guías de ondas para otros instrumentos reales, explicaremos la secuenciación de sonido para dar cuerpo a la síntesis de sonidos y poder explotarla como fuente musical, y veremos métodos de procesado de sonido, aplicables tanto a los sonidos que estamos creando, como a cualquier grabación de audio (eco, reverberación, delay, compresión, etc).
