scm_core.asv

3GPP推荐的SCM信道模型的MATLAB源码

%SCM_CORE Channel coefficient computation for a geometric channel model
%   [H DELTA_T FINAL_PHASES FINAL_PHASES_LOS]=SCM_CORE(SCMPAR, LINKPAR,
%   ANTPAR, BULKPAR, BSGAIN, BSGAIN_LOS, MSGAIN, MSGAIN_LOS, OFFSET_TIME,
%   BSGAINISSCALAR, MSGAINISSCALAR) This is the scm_core aka the big for
%   loop. It implements the formulas in [1, Sec. 5.4] and [1, Sec. 5.5].
%
%   Outputs:
%
%   H               - [UxSxNxTxK] array of channel coefficients
%   DELTA_T         - time sampling intervals (in seconds) for all links
%   FINAL_PHASES    - final phases of all subpaths in degrees over (-180,180)
%   FINAL_PHASES_LOS- final phases for LOS paths in degrees over (-180,180)
%   
%   Inputs:
%
%   SCMPAR          - input struct, see SCMPARSET
%   LINKPAR         - input struct, see LINKPARSET
%   ANTPAR          - input struct, see ANTPARSET
%   BULKPAR         - input BULKPAR, see GENERATE_BULK_PAR
%   BSGAIN          - [KxSxNxM] array of interpolated antenna field
%                     patterns (complex)
%   BSGAIN_LOS      - [KxS] array of interpolated antenna field patterns
%                     (complex) for LOS paths. Only used with the LOS
%                     option; it is set to scalar otherwise.
%   MSGAIN          - [KxUxNxM] array of interpolated antenna field
%                     patterns (complex)
%   MSGAIN_LOS      - [KxU] array of interpolated antenna field patterns
%                     (complex) for LOS paths. Only used with the LOS
%                     option; it is set to scalar otherwise.
%   OFFSET_TIME     - time offset added to the initial phase (set to zero by default)
%   BSGAINISSCALAR  - this is 1 if BsGain is uniform over azimuth, 0 otherwise.
%   MSGAINISSCALAR  - this is 1 if MsGain is uniform over azimuth, 0 otherwise.
%
%   With 'polarized' option:
%
%   BSGAIN          - [KxSx2xNxM] array of interpolated antenna field
%                     patterns (complex), where the third dimension are
%                     the patterns for [V H] polarizations. 
%   MSGAIN          - [KxUx2xNxM] array of interpolated antenna field
%                     patterns (complex), where the third dimension are
%                     the patterns for [V H] polarizations. 
%
%   To compile the ANSI-C written optimized core, type
%
%       mex scm_mex_core.c
%
%   at MATLAB prompt. For further documentation on the ANSI-C implementation 
%   of the SCM_CORE, see SCM_MEX_CORE.

%   Authors: Giovanni Del Galdo (TUI), Jussi Salmi (HUT), Marko Milojevic (TUI), 
%   Christian Schneider (TUI), Jari Salo (HUT), Pekka Ky鰏ti (EBIT), 
%   Daniela Laselva (EBIT)
%   $Revision: 0.3 $  $Date: Jan 13, 2005$


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%      %%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%            %%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%                  %%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%                        %%%%%%%%%%%%%%
%%%%%%%%%%%                              %%%%%%%%%%%
%%%%%%%%                                    %%%%%%%%
%%%%%                                          %%%%%
%%                   --------                     %%
function [H, delta_t, output_SubPathPhases, output_Phi_LOS] = scm_core (scmpar,linkpar,antpar,bulkpar,BsGain,BsGain_Theta_BS,MsGain,MsGain_Theta_MS,offset_time, BsGainIsScalar, MsGainIsScalar)
%%                   --------                     %%
%%%%%                                          %%%%%
%%%%%%%%                                    %%%%%%%%
%%%%%%%%%%%                              %%%%%%%%%%%
%%%%%%%%%%%%%%                        %%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%                  %%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%            %%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%      %%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% scm_core AKA the big for loop. It implements the formula in 5.4
%
%
% offset_time [samples] = defines the starting point (in samples) of the
%                         time axis
% Examples: you want to calculate 1000 time samples calling the scm_core
%           twice (everytime for 500 timesamples)
%           H1 = scm_core (..., 0)
%           H2 = scm_core (..., 500)
%
%   Revision history:
%  
% ______09-June-2004_____
% 1- now scm_core uses path_powers and not subpath_powers
% 2- fixed PHI_LOS
% 3- fixed the antenna gains for the LOS component
%
% ______15-June-2004_____
% 
% 1- the powers must be in linear scale
% 2- LOS probability has been removed since it's never used
% 3- the shadow fading will be applied outside the scm_core, thus it's
%    removed
%
% ______16-June-2004_____
% 
% 1- removed the string argument to call 'safe' 'sparse' etc.  Now it
%    always runs the 'safe' way
% 2- Giovanni: added the polarized option
% 3- removed the sqrt from the gains - they are complex and expressed in
%    tension
% 4- added the phase output
%
% ______21-June-2004_____
% 1- fixed a bug in the polarized loop. du and aoas were missing
%
% ______21-June-2004_____
% 1- fixed a bug in the polarized loop. aoas and aods are degrees! 
%
% ______19-July-2004_____
% 1- Jari: help text and minor editorial work
%
% ______26-July-2004_____ (Jari)
% 1- changed 'polarized' option so that mean power is normalized
% 2- changed the way S and U are determined
% 3- changed MsElementSpacingULA to MsElementPosition
% 4- added ANSI-C core
%
% ______29-July-2004_____ (Jari)
% 1- Added two new input args: BsGainIsScalar and MsGainIsScalar
%
% ______14-September-2004_____ (Jari)
% 1- Modifications for intra-cluster delay spread (ANSI-C part only)
% 2- Fixed delta_t so that SampleDensity means samples per ? wavelength
%
%______21-Sep-2004_____ (Giovanni)
% 1- adding the intra-cluster DS for all cases
%
%______7-Dec-2004_____ (Jari)
% 1- Bug fix in LOS option: MsElementPosition(s) changed to
% BsElementPosition(s) in a for loop
% 2- some modifications for the public release
%
%______14-Jan-2005_____ 
% 1- (Marko) Added the option XpdIndependentPower
% 2- (Jussi) Added the option XpdIndependentPower for ANSI-C core
%


DEBUG_MODE_FLAG   = 1;
PROFILE_MODE_FLAG = 0;
DISPLAY_MODE_FLAG = 0;

S=size(BsGain,2);               % number of receiving antennas
U    = size(MsGain,2);               % number of transmitting antennas
N    = scmpar.NumPaths;           % number of paths
T    = scmpar.NumTimeSamples;     % number of time samples
K    = length(linkpar.MsNumber);  % number of links
M    = scmpar.NumSubPathsPerPath; % number of subpaths


scmpar.IntraClusterDsUsed='no';    % this is fixed in this version
if strcmpi(scmpar.IntraClusterDsUsed,'yes')
    LM = bulkpar.MidPathOrder;
    LN = bulkpar.NumSubPathsPerMidpath;
    L = length(LN);
    
    path_powers_all = bulkpar.path_powers_all;
    
    if(sum(LN) ~= M)
        error('Sum over NumSubPathsPerMidpath must equal NumSubPathsPerPath.');
    end
else    % special case, one midpath only
    LM = 1:M;
    LN=M;
    L = length(LN);
    
    path_powers_all = bulkpar.path_powers;
end


H = zeros(U,S,N,T,K);


% define element spacing vectors if scalars are given
if (length(antpar.MsElementPosition)==1)
    antpar.MsElementPosition=[0:antpar.MsElementPosition:antpar.MsElementPosition*(U-1)];
end

if (length(antpar.BsElementPosition)==1)
    antpar.BsElementPosition=[0:antpar.BsElementPosition:antpar.BsElementPosition*(S-1)];
end


if DISPLAY_MODE_FLAG
    disp('  ');disp('  ');disp('  ')
    disp('          _______________________ ');
    disp('         ''                       ''');
    disp('         |                       |');
    disp('         |  Welcome              |');
    disp('         |      to the SCM core  |');
    disp('         |                       |');
    disp('         |                       |');    
    disp('         ''-----------------------''');    
    disp('  ')
    disp('  ')
    disp('  ')
end    


% Set internal parameters
speed_of_light=2.99792458e8;
wavelength=speed_of_light/scmpar.CenterFrequency;

% dummy
output_Phi_LOS       = zeros(K,1);

% let's make the time axis - for that we need to check UniformTimeSampling
% and the MSs' velocities
% Note: SampleDensity is samples per half wavelength.
if strcmp(scmpar.UniformTimeSampling,'yes')
    
    max_vel = max(linkpar.MsVelocity);
    delta_t = repmat((wavelength / max_vel)/2/scmpar.SampleDensity,K,1);
    
else % 'UniformTimeSampling' is 'no'
    
    delta_t = (wavelength ./ linkpar.MsVelocity.')./2/scmpar.SampleDensity ;
    
end
t = repmat(delta_t,1,T).*repmat([0:T-1]+offset_time,K,1); % matrix containing the time axes for all links [KxT]
% t = repmat(delta_t,1,T).*repmat([1:T],K,1); % matrix containing the time axes for all links [KxT]

%%%%%%
if DEBUG_MODE_FLAG
    figure
    plot(t.')
    xlabel('samples')
    ylabel('time [sec]')
    grid on
    delta_t
end
%%%%%%

%%%%%%
if PROFILE_MODE_FLAG
    profile on -detail builtin -history
end
%%%%%%



k_CONST = 2*pi/wavelength;

%%%%%%%%%%%%%%%%%%ANSI-C core part%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%check if ANSI-C core is used
if strcmp(lower(scmpar.AnsiC_core),'yes')
    
    % different modes
    GENERAL = 1;
    POLARIZED = 2;
    LOS = 3;   
    
    look_up_points = scmpar.LookUpTable; % set if look-up table is used for sin/cos
    
    
    if (BsGainIsScalar && MsGainIsScalar)
        GainsAreScalar = 1;
    else 
        GainsAreScalar = 0;
    end
    
    if ~strcmpi(scmpar.ScmOptions,'polarized')
        
        if DISPLAY_MODE_FLAG
            disp('entering main loop...');
        end
        
        % adjusting parameters for calling the C routine
        d_u = antpar.MsElementPosition*wavelength;
        d_s = antpar.BsElementPosition*wavelength;
        aod = bulkpar.aods(1:K,1:N,1:M)*pi/180;
        aoa = bulkpar.aoas(1:K,1:N,1:M)*pi/180;  
        phase = bulkpar.subpath_phases(1:K,1:N,1:M)*pi/180;
        v = linkpar.MsVelocity(1:K);
        theta_v = linkpar.MsDirection(1:K)*pi/180;
        sq_Pn = sqrt(path_powers_all(1:K,1:N*L));
        
        % calling the C-mex routine for general coefficients 
        [H output_SubPathPhases] = scm_mex_core(GENERAL, BsGain, MsGain, aod, aoa, d_s, d_u, phase, t, k_CONST, v, theta_v, sq_Pn, look_up_points, U, S, N, M, K, T, GainsAreScalar);
        
        %output_SubPathPhases = prin_value((output_SubPathPhases*180/pi + bulkpar.subpath_phases)); %changed due tests
        output_SubPathPhases = prin_value((output_SubPathPhases*180/pi + bulkpar.subpath_phases(1:K,1:N,1:M)));
        
    else % it's polarized!
        if DISPLAY_MODE_FLAG
            disp('entering polarized option...');
        end    
        
        output_SubPathPhases = zeros(K,4,N,M);
        %temp_output_SubPathPhases = zeros(K,N,M);
        
        % BsGain must have size: [K x S x 2 x N x M]
        % the first dimension in the polarization must be vertical
        %
        % MsGain must have size: [K x U x 2 x N x M]
        % the first dimension in the polarization must be vertical
        %
        % subpath_phases has size: [K x 4 x N x M]
        % check that vv and vh etc are the ones we think they are 
        %
        % bulkpar.xpd has size: [K,2,N]
        % keyboard
        % temp = zeros(M,T);    
        %     keyboard
        % adjusting parameters for calling the C routine
        d_u = antpar.MsElementPosition * wavelength;