循环功率谱 Matlab Simulink仿真程序

#define S_FUNCTION_NAME Channel #define S_FUNCTION_LEVEL 2 #include "simstruc.h" #include <math.h> #include <stdio.h> enum {ChnlType = 0, SampleFreq, Fc, Seed, ChnlIndex, NUM_ARGS}; #define CHNL_TYPE_ARG(S) (ssGetSFcnParam(S, ChnlType)) #define SAMPLE_FREQ_ARG(S) (ssGetSFcnParam(S, SampleFreq)) #define FC_ARG(S) (ssGetSFcnParam(S, Fc)) #define SEED_ARG(S) (ssGetSFcnParam(S, Seed)) #define CHNL_IDX(S) (ssGetSFcnParam(S, ChnlIndex)) #define GET_CHNL_TYPE(S) (mxGetPr(CHNL_TYPE_ARG(S)))[0] #define GET_SAMPLE_FREQ(S) (mxGetPr(SAMPLE_FREQ_ARG(S)))[0] #define GET_FC(S) (mxGetPr(FC_ARG(S)))[0] #define GET_SEED(S) (mxGetPr(SEED_ARG(S)))[0] #define GET_CHNL_IDX(S) (mxGetPr(CHNL_IDX(S)))[0] #define PI 3.1415926536 #define SQRT2 1.4142135623 #define MULTI 16807 // 7^5 #define MODE 2147483647 // 2^31-1 #define MAX_ENERGY 1e10 #define MAX_CHNL_SMPL 150 #define MAX_CHNL_LEN 250 // Unit: ns #define TOTAL_CHNL_NUM 100 #define SEL_CHNL_NUM 90 #define FFT_SIZE 64 #define CHNL_OUT_SIZE 32 // 参考输出信道长度 #define WIN_WIDTH_CM1 6 #define WIN_WIDTH_CM2 8 #define WIN_WIDTH_CM3 12 #define WIN_WIDTH_CM4 16 // real_T chnlR[SEL_CHNL_NUM][MAX_CHNL_SMPL], chnlI[SEL_CHNL_NUM][MAX_CHNL_SMPL]; // int_T chnlLen[SEL_CHNL_NUM]; // static int_T chnlIndex; int_T FindStart(real_T* dataR, real_T *dataI, int_T chnlLen, SimStruct *S) { int_T ii, jj; real_T maxEnergy = 0; real_T energy = 0; int_T start = 0; int_T chnlType=(int_T)GET_CHNL_TYPE(S); int_T winWidth; switch(chnlType) { case 1: winWidth = WIN_WIDTH_CM1; break; case 2: winWidth = WIN_WIDTH_CM2; break; case 3: winWidth = WIN_WIDTH_CM3; break; case 4: winWidth = WIN_WIDTH_CM4; break; default: winWidth = WIN_WIDTH_CM3; } for(ii=0;ii<chnlLen-winWidth;ii++) { energy = 0; for (jj=0;jj<winWidth;jj++) { energy += dataR[ii+jj] * dataR[ii+jj] + dataI[ii+jj] * dataI[ii+jj]; } if (energy>maxEnergy) { start = ii; maxEnergy = energy; } } return start; } real_T Rand01(int_T seed, SimStruct *S) { static real_T x; real_T y,t,n; int_T calledFlag = ssGetIWork(S)[0]; // load calledFlag if (calledFlag == 0) { ssSetIWorkValue(S, 0, 1); // update calledFlag x = seed; } t = x*MULTI/MODE; y = modf(t,&n); x = x*MULTI - n*MODE; return x/MODE; } /***************************************************** real_T RandGaussian(int_T seed) 函数功能：产生一个高斯分布的随机变量，其均值为0，方差为1 返回值：产生的高斯随机变量，real_T型 参数说明： seed: 产生随机数的种子 其它说明 修改日期 版本号 修改人 修改内容 ------------------------------------------------------ 2004/07/13 V1.0 张在琛 创建 ******************************************************/ real_T RandGaussian(int_T seed, SimStruct *S) { real_T rand01X, rand01Y; // 0 - 1之间均匀分布的随机数 real_T randG; // 高斯分布的随机数 rand01X = Rand01(seed, S); rand01Y = Rand01(seed, S); randG = sqrt(-2*log(rand01X))*cos(2*PI*rand01Y); return randG; } static void checkPara(SimStruct *S) { real_T chnlType, sampleFreq, fc; chnlType = (int_T)GET_CHNL_TYPE(S); sampleFreq = GET_SAMPLE_FREQ(S); // Unit: MHz fc = GET_FC(S); // Unit: MHz if ((chnlType < 0)||(chnlType > 4)) { ssSetErrorStatus(S,"Parameter CHNL_TYPE error!\n"); return; } if (1000 / sampleFreq * MAX_CHNL_SMPL < MAX_CHNL_LEN) { ssSetErrorStatus(S,"Warning: MAX_CHNL_SMPL may be too small for CM4!\n"); return; } if (fc < 0) { ssSetErrorStatus(S,"Parameter FC error!\n"); return; } } /***************************************************** void GetThEnergy(real_T *pThEnergy, real_T *pChnlEnergy, SimStruct *S) 函数功能：找到选择信道的门限能量值 返回值：void型： 参数说明： pThEnergy: 找到的能量门限值 pChnlEnergy: 每一个信道冲激响应序列的能量 其它说明：在产生的TOTAL_CHNL_NUM个信道中，需要选取SEL_CHNL_NUM个能量 最大的信道最为最后的输出。本函数寻找能量的门限值，能量大于 等于此门限值的信道为最后输出的信道。因为能量值为real_T型， 我们认为两个信道能量相等的情况极少出现。为简单起见，不考虑 两个信道能量相等的情况。如果出现了这一情况，则channel模块 中可能会出现选取了能量较小的信道的情况。即时出现了这一情况， 也只是对性能有细微的影响，不会影响程序的正常运行。 ******************************************************/ void GetThEnergy(real_T *pThEnergy, real_T *pAvgE, real_T *pChnlEnergy, SimStruct *S) { int_T ii, jj; // 循环变量 int_T index; // 指示最小能量值位置的变量 int_T rmvNum; // 需要去除的能量较小的信道的数目 int_T pFlag[TOTAL_CHNL_NUM]; // 标志数组，指示对应位置的信道能量是否 // 已被选为过最小值，如是，设为1，否则为0 real_T miniEnergy; // 最小的能量值 real_T tempE; rmvNum = TOTAL_CHNL_NUM - SEL_CHNL_NUM; // 检查常量 if (rmvNum < 0) { ssSetErrorStatus(S,"GetThEnergy detects TOTAL_CHNL_NUM < SEL_CHNL_NUM!\n"); return; } for (ii=0; ii<TOTAL_CHNL_NUM; ii++) { if (pChnlEnergy[ii] > MAX_ENERGY) { ssSetErrorStatus(S,"GetThEnergy detects MAX_ENERGY is too small!\n"); return; } } // 初始化标志数组 for (ii=0; ii<TOTAL_CHNL_NUM; ii++) { pFlag[ii] = 0; } // 循环(rmvNum+1)，依次找出(rmvNum+1)最小的能量值，最后一个即为输出 for (ii=0; ii<=rmvNum; ii++) { miniEnergy = MAX_ENERGY; for (jj=0; jj<TOTAL_CHNL_NUM; jj++) { if ((pChnlEnergy[jj] < miniEnergy) && (pFlag[jj] == 0)) { miniEnergy = pChnlEnergy[jj]; index = jj; } } pFlag[index] = 1; // 标志设为1，避免下次再被选为最小 } // 以下为计算选择后信道的归一化因子 tempE = 0; for (ii=0; ii<TOTAL_CHNL_NUM; ii++) { if (pChnlEnergy[ii]>=miniEnergy) { tempE += pChnlEnergy[ii]; } } //*pAvgE = tempE / SEL_CHNL_NUM; *pAvgE = 1; *pThEnergy = miniEnergy; } /***************************************************** void ChnlInit(int_T chnlType, real_T sampleFreq, real_T fc, int_T seed, SimStruct *S) 函数功能：产生一组信道冲激响应 返回值：void型： 参数说明： chnlType: 信道的类型（有效值为1 - 4） sampleFreq: 信道的采样频率（单位：MHz） fc: 载波中心频率（单位：MHz）（3696 MHz或4224 MHz） seed: 随机数的种子 其它说明 ******************************************************/ void ChnlInit(int_T chnlType, real_T sampleFreq, real_T fc, int_T seed, SimStruct *S) { real_T lambdaCluster; // 簇的到达率（单位：GHz） real_T lambdaRay; // 径的到达率（单位：GHz） real_T gammaCluster; // 簇的衰减因子（单位：ns） real_T gammaRay; // 径的衰减因子（单位：ns） real_T stdLn1; // 簇衰减中对数正态变量的标准差 real_T stdLn2; // 径衰减中对数正态变量的标准差 int_T nLos; // 0表示视线（LOS）传输，1表示非视线传输（NLOS） real_T stdShdw; // 整个冲激响应的对数正态衰落的标准差 real_T stdLamCluster; // 簇到达间隔的标准差（单位：ns） real_T stdLamRay; // 径到达间隔的标准差（单位：ns） real_T muConst; // 中间变量 real_T lnXi; real_T tVal; // 时间值 real_T hVal; // 信道响应值 real_T mu; real_T lnBeta; int_T pk; // pk = 1或-1，代表信号的反射 int_T ii, kk, k0; // 循环变量 real_T temp, temp1, temp2; // 临时变量 real_T chnlTmp[TOTAL_CHNL_NUM][MAX_CHNL_SMPL]; int_T chnlLenTmp[TOTAL_CHNL_NUM]; real_T chnlEnergy[TOTAL_CHNL_NUM];