%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%% Function CIRpowers.m                                         %%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% [tap_variances,L]=CIRpowers(Model)
%% This function returns the taps variances of a static multipath channel
%% and the channel order

%% implements sample-spaced channel filtering, no signal interpolation.

%% delays : length L vector of tap delays (in microseconds) (L: number of taps)
%% powers : length L vector of tap powers (in dB) (L: number of taps)
%%          (note that the coefficients will automatically normalized so that
%%          the channel has an overall gain of 0 dB)
%%
%% Octobre 2000. Bertrand Muquet, Sebastien Simoens, Shengli Zhou

function [tap_variances,L]=CIRpowers(Model)

  T=5e-2;

  switch Model
    case 1
      delays=1e-3*[0 10 20 30 40 50 60 70 80 90 110 140 170 200 240 290 340 390];powers=[0 -0.9 -1.7 -2.6 -3.5 -4.3 -5.2 -6.1 -6.9 -7.8 -4.7 -7.3 -9.9 -12.5 -13.7 -18.0 -22.4 -26.7];
    case 2
      delays=1e-3*[0 10 20 30 50 80 110 140 180 230 280 330 380 430 490 560 640 730];powers=[-2.6 -3.0 -3.5 -3.9 0.0 -1.3 -2.6 -3.9 -3.4 -5.6 -7.7 -9.9 -12.1 -14.3 -15.4 -18.4 -20.7 -24.6];
    case 3
      delays=1e-3*[0 10 20 30 50 80 110 140 180 230 280 330 400 490 600 730 880 1050];powers=[-3.3 -3.6 -3.9 -4.2 0.0 -0.9 -1.7 -2.6 -1.5 -3.0 -4.4 -5.9 -5.3 -7.9 -9.4 -13.2 -16.3 -21.2];
    case 5
  delays=1e-3*[0 10 20 40 70 100 140 190 240 320 430 560 710 880 1070 1280 1510 1760];powers=[-4.9 -5.1 -5.2 -0.8 -1.3 -1.9 -0.3 -1.2 -2.1 0.0 -1.9 -2.8 -5.4 -7.3 -10.6 -13.4 -17.4 -20.9];
  end

  %%The delays are normalized
  sz=size(delays);
  if (and(sz(1) ~= 1,sz(2) == 1)) delays=delays.';
  elseif (and(sz(1) ~= 1,sz(2) ~= 1)) 'Error: delays must be a vector';
  end
  delays=delays/T; %Delays are expressed in number of samples.
  
  nbtaps=length(powers);
  len_cir=1+round(max(delays));
  tap_variances=zeros(1,len_cir);
  
  %%Compute the amplitude of each tap
  sz=size(powers);
  if (and(sz(1) ~= 1,sz(2) == 1)) powers=powers.';
  elseif (and(sz(1) ~= 1,sz(2) ~= 1)) 'Error: powers must be a vector';
  end

  %% Powers are in dB -> computes their variances to recombine power taps.
  variances=10.^(powers/10);
  %% powers are normalized
  variances=variances/sum(variances);

  %%Then the discrete-time CIR is computed by rounding each tap
  %%to the next sample
  for i=1:nbtaps
    tap_variances(1+round(delays(i)))=tap_variances(1+round(delays(i)))+ variances(i);
  end

  L=length(tap_variances)-1;
  

