How can I transform a time signal into a third octave spectrum? I try to make a third-octave spectrum from a time signal and I have some troubles with the algorithm. I'm not an expert in signal processing so I don't know what is wrong.

Here what my algorithm looks like :

Time signal -> Z transform of the signal -> First octave band filtering -> Quadratic sum of the result's values => First octave band value.
Under factor 2 sampling ->Z transform of the signal -> Second octave band filtering -> Quadratic sum of the result's values =>Second octave band.

Now I have this, I try to test it on Matlab with a time signal about 2 millions points. 

I have this code to calculate all the octave's bands :

    fc = fMin*2.^((0:bmax)*bw);             %Center frequencies.
    fl = fc*2^(-bw/2);                      %Lower cutoffs.
    fu = fc*2^(+bw/2);                      %Upper cutoffs.


And this one to calculate the Z transform of the time signal :

    x = [1 2 3];
    syms z;
    d = z.^(-1*(0:numel(x)-1));
    zTransform = sum(x.*d)


My main problem is, I'm not completely sure that my algorithm is right and if it is, how can I filter the z transform equation that the code above gives me with each third-octave band filters?

Thanks in advance, it's been a while I'm on this problem.
maybe either the bandpass filters or the peaking EQ filters in the audio EQ cookbook can help.  bandwidth is defined in octaves in that.  the output, $y$, of each 1/3 octave BPF filter should be squared and low-pass filtered to get the mean square, $\overline{y^2}$.  and then 

$$dB = 10 \log_{10}\Big(\,\overline{y^2}\,\Big)$$

is applied to that to get dB.

are you planning on doing this to the whole audible spectrum?  31 bands?