信道容量的迭代算法实现

   1: #include 
   
   

   2: #include 
   
   

   3: #include 
   
   

   4: #include 
   
   

   5: using namespace std;

   6:  

   7: #define FLOAT_MINUS_PRECISION 0.00001

   8: typedef vector<float*> VEC_PFLOAT;

   9:  

  10: //迭代计算信道容量,参数值为信源,信宿符号个数和信道转移概率矩阵,返回信道容量

  11: float GetCapacity(int nSourceSymbol,int nHostSymbol,const VEC_PFLOAT& vTransMatrix)

  12: {

  13:     //信道容量初始化为最小值

  14:     float fCapacity = FLT_MIN;

  15:     //信源概率分布

  16:     float *pfSoureProb = new float[nSourceSymbol];

  17:     //初始化信源分布为均匀分布

  18:     int i;

  19:     for (i = 0; i < nSourceSymbol; i++)

  20:     {

  21:         pfSoureProb[i] = 1.0 / nSourceSymbol;

  22:     }

  23:     //初始化φ函数

  24:     VEC_PFLOAT vPhi;

  25:     for (i = 0; i < nSourceSymbol; i++)

  26:     {

  27:         float *pfTemp = new float[nHostSymbol];        

  28:         vPhi.push_back(pfTemp);

  29:     }

  30:     //设置精度;

  31:     const float cfDelta = 0.02f;

  32:     float fPrecision;

  33:     //迭代计算

  34:     int j,k;

  35:     float *pfSum = new float[nSourceSymbol];

  36:     do 

  37:     {

  38:         for (i = 0; i < nSourceSymbol; i++)

  39:         {

  40:             for (j = 0; j < nHostSymbol; j++)

  41:             {                

  42:                 //计算ΣPi*Pji

  43:                 float fSum = 0.0f;

  44:                 for (k = 0; k < nSourceSymbol; k++)

  45:                 {

  46:                     fSum += pfSoureProb[i] * vTransMatrix[k][i];

  47:                 }

  48:                 vPhi[i][j] = pfSoureProb[i] * vTransMatrix[j][i] / fSum;

  49:             }            

  50:         }

  51:         float fSumDeno = 0.0f;    //分母求和                

  52:         for (i = 0; i < nSourceSymbol; i++)

  53:         {        

  54:             float fSum = 0.0f;

  55:             for (j = 0; j < nHostSymbol; j++)

  56:             {

  57:                 fSum += vTransMatrix[j][i] * logf(vPhi[i][j]);

  58:             }

  59:             pfSum[i] = expf(fSum);

  60:             fSumDeno += pfSum[i];

  61:         }

  62:         for (i = 0; i < nSourceSymbol; i++)

  63:         {            

  64:             pfSoureProb[i] = pfSum[i] / fSumDeno;

  65:         }

  66:         //计算新一轮的容量

  67:         float fNewC = logf(fSumDeno);

  68:         //计算精度

  69:         fPrecision = fabs(fNewC - fCapacity) / fCapacity;

  70:         fCapacity = fNewC;

  71:     } while(fPrecision - cfDelta > 0.0f);

  72:     //释放临时资源...

  73:     delete []pfSum;

  74:     for (i = 0; i < vPhi.size(); i++)

  75:     {

  76:         float* pfTemp = vPhi.at(i);

  77:         delete pfTemp;

  78:     }

  79:     vPhi.clear();

  80:     return fCapacity;

  81: }

  82: int main()

  83: {

  84:     //转移矩阵

  85:     VEC_PFLOAT vTransMatrix;

  86:     int nCol,nLine;

  87:     cout<<"请输入信源符号个数:";

  88:     cin>>nLine;

  89:     cout<<"请输入信宿符号个数:";

  90:     cin>>nCol;

  91:     cout<<"请依次输入"<
   
   
    
    <<
    
    "行信道转移概率矩阵:(以空格隔开每个概率)/n";
   
   

  92:     for (int i = 0; i < nLine; i++)

  93:     {

  94:         float *pfTemp = new float[nCol];

  95: Label1:

  96:         float fSum = 0.0f;    

  97:         cout<<"X"<<<":";

  98:         for (int j = 0; j < nCol - 1; j++)

  99:         {

 100:             cin>>pfTemp[j];

 101:             fSum += pfTemp[j];

 102:         }

 103:         if (1.0f - fSum < 0)

 104:         {

 105:             cout<<"转移概率和应该为1,请重新输入!/n";

 106:             goto Label1;

 107:         }

 108:         else

 109:         {

 110:             pfTemp[j] = 1.0f - fSum;

 111:             cout<<"信源符号X"<<<"的转移概率分别为:";

 112:             for(int k = 0; k < nCol; k++)

 113:                 cout<

 114:             cout<

 115:         }

 116:         vTransMatrix.push_back(pfTemp);

 117:     }

 118:     cout<<"信道容量为:"<

 119:     for (int k = 0; k < vTransMatrix.size(); k++)

 120:     {

 121:         float* pfTemp = vTransMatrix.at(k);

 122:         delete pfTemp;

 123:     }

 124:     vTransMatrix.clear();

 125:     return 0;

 126: }

算法描述请参照信息论相关理论知识