SoloCodigo

Hola espero que esta respuesta te ayude.
En primer lugar decirte que la fft es una función de matlab que puede que desaparezca, pues existen funciones mejores para obtener la densidad espectral de potencia de una señal.

Lo primero que tienes que entender es qué es tu función wav: cuando haces un wavread

[cancion, fs, nbits]=wavread('MiPCMis Documentoscancion.wav');

Lo que obtienes es simplemente una serie de muestras de la canción (cancion), también conoces a qué frecuencia de muestreo (fs) se ha tomado cada muestra de la canción cuando se creó ese wav y el número de bits (nbits) que se usan para codificar cada muestra.

La función cancionfft=fft(canción,nfft) devuelve la transformada de fourier de la señal canción calculada en tantas frecuencias como le indiques con el número nfft, es decir, tú puedes ajustar la precisión con la que va a calcular la transformada, evidentemente hasta el límite que le permite la información que tiene. La función fft tiene para un principiante algunos problemas: el primero es que si haces un plot(cancionfft) verás que obtienes algo simétrico, es decir, ha calculado la parte positiva y la negativa del espectro. Sabes que el centro de simetría es la frecuencia cero, pero estás tratando de analizar un wav y la parte negativa del espectro no te interesa para nada. Por otra parte, para saber a qué frecuencia corresponde cada muestra calculada has de tener en cuenta la nfft y la fs y calcularla.  
Por último la transformada no es la densidad de potencia de la señal, en realidad, como un wav no tiene ninguna referencia de potencia real, calcular la potencia exacta pero si puedes toda la información espectral relativa. Me estoy liando.

Existe otro camino más sencillo que te recomiendo para obtener lo que necesitas, y es utilizar la función Pwelch, que calcula directamente la densidad espectral de potencia de tu señal wav. Pwelch hace una fft modificada, no es difícil pero es algo largo de explicar. Te permite ajustar bastantes parámetros, y algunos los ajusta por defecto, y obtiene directamente la densidad espectral de potencia de la señal (la transformada al cuadrado). Los siguientes pasos te servirán:

[cancion, fs, nbits]=wavread('MiPCMis Documentoscancion.wav');
[DEPcancion,eje_frecuencias]=Pwelch(cancion,[],[],[],fs);

ahora te recomiendo que lo pases a decibelios para verlo más claro, si supieras que las unidades de potencia de tus muestras wav codifican unidades de Pascales tendrías que hacer

DEPcancion_dB=10*log10(DEPcancion./4e-10);

En realidad no lo sabes, pero si lo haces (aunque no dividas por la referencia 4e-10 eso da lo mismo), verás el espectro mucho más claro.

Ahora haz
plot(eje_frecuencias,DEPcancion_dB);
ylabel('dB');
xlabel('Hercios');
title('DEP cancion');

Y ya está, el eje x es el eje de frecuencia real en Hercios correspondiente a cada muestra, y en el eje tienes la forma de la DEP en unidades logarítmicas.

Espero que te saque del apuro, mucha suerte!