git.gag.com Git - debian/gnuradio/blob - gr-trellis/src/lib/fsm.cc

   1 /* -*- c++ -*- */
   2 /*
   3  * Copyright 2002 Free Software Foundation, Inc.
   4  *
   5  * This file is part of GNU Radio
   6  *
   7  * GNU Radio is free software; you can redistribute it and/or modify
   8  * it under the terms of the GNU General Public License as published by
   9  * the Free Software Foundation; either version 2, or (at your option)
  10  * any later version.
  11  *
  12  * GNU Radio is distributed in the hope that it will be useful,
  13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15  * GNU General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU General Public License
  18  * along with GNU Radio; see the file COPYING.  If not, write to
  19  * the Free Software Foundation, Inc., 51 Franklin Street,
  20  * Boston, MA 02110-1301, USA.
  21  */
  22
  23 #include <cstdio>
  24 #include <string>
  25 #include <iostream>
  26 #include <fstream>
  27 #include <stdexcept>
  28 #include <cmath>
  29 #include "base.h"
  30 #include "fsm.h"
  31
  32
  33 fsm::fsm()
  34 {
  35   d_I=0;
  36   d_S=0;
  37   d_O=0;
  38   d_NS.resize(0);
  39   d_OS.resize(0);
  40   d_PS.resize(0);
  41   d_PI.resize(0);
  42   d_TMi.resize(0);
  43   d_TMl.resize(0);
  44 }
  45
  46 fsm::fsm(const fsm &FSM)
  47 {
  48   d_I=FSM.I();
  49   d_S=FSM.S();
  50   d_O=FSM.O();
  51   d_NS=FSM.NS();
  52   d_OS=FSM.OS();
  53   d_PS=FSM.PS(); // is this going to make a deep copy?
  54   d_PI=FSM.PI();
  55   d_TMi=FSM.TMi();
  56   d_TMl=FSM.TMl();
  57 }
  58
  59 fsm::fsm(int I, int S, int O, const std::vector<int> &NS, const std::vector<int> &OS)
  60 {
  61   d_I=I;
  62   d_S=S;
  63   d_O=O;
  64   d_NS=NS;
  65   d_OS=OS;
  66
  67   generate_PS_PI();
  68   generate_TM();
  69 }
  70
  71 //######################################################################
  72 //# Read an FSM specification from a file.
  73 //# Format (hopefully will become more flexible in the future...):
  74 //# I S O (in the first line)
  75 //# blank line
  76 //# Next state matrix (S lines, each with I integers separated by spaces)
  77 //# blank line
  78 //# output symbol matrix (S lines, each with I integers separated by spaces)
  79 //# optional comments
  80 //######################################################################
  81 fsm::fsm(const char *name)
  82 {
  83   FILE *fsmfile;
  84
  85   if((fsmfile=fopen(name,"r"))==NULL)
  86     throw std::runtime_error ("fsm::fsm(const char *name): file open error\n");
  87     //printf("file open error in fsm()\n");
  88
  89   fscanf(fsmfile,"%d %d %d\n",&d_I,&d_S,&d_O);
  90   d_NS.resize(d_I*d_S);
  91   d_OS.resize(d_I*d_S);
  92
  93   for(int i=0;i<d_S;i++) {
  94     for(int j=0;j<d_I;j++) fscanf(fsmfile,"%d",&(d_NS[i*d_I+j]));
  95   }
  96   for(int i=0;i<d_S;i++) {
  97     for(int j=0;j<d_I;j++) fscanf(fsmfile,"%d",&(d_OS[i*d_I+j]));
  98   }
  99
 100   generate_PS_PI();
 101   generate_TM();
 102 }
 103
 104
 105
 106
 107 //######################################################################
 108 //# Automatically generate the FSM from the generator matrix
 109 //# of a (n,k) binary convolutional code
 110 //######################################################################
 111 fsm::fsm(int k, int n, const std::vector<int> &G)
 112 {
 113
 114   // calculate maximum memory requirements for each input stream
 115   std::vector<int> max_mem_x(k,-1);
 116   int max_mem = -1;
 117   for(int i=0;i<k;i++) {
 118     for(int j=0;j<n;j++) {
 119       int mem = -1;
 120       if(G[i*n+j]!=0)
 121         mem=(int)(log(G[i*n+j])/log(2.0));
 122       if(mem>max_mem_x[i])
 123         max_mem_x[i]=mem;
 124       if(mem>max_mem)
 125         max_mem=mem;
 126     }
 127   }
 128
 129 //printf("max_mem_x\n");
 130 //for(int j=0;j<max_mem_x.size();j++) printf("%d ",max_mem_x[j]); printf("\n");
 131
 132   // calculate total memory requirements to set S
 133   int sum_max_mem = 0;
 134   for(int i=0;i<k;i++)
 135     sum_max_mem += max_mem_x[i];
 136
 137 //printf("sum_max_mem = %d\n",sum_max_mem);
 138
 139   d_I=1<<k;
 140   d_S=1<<sum_max_mem;
 141   d_O=1<<n;
 142
 143   // binary representation of the G matrix
 144   std::vector<std::vector<int> > Gb(k*n);
 145   for(int j=0;j<k*n;j++) {
 146     Gb[j].resize(max_mem+1);
 147     dec2base(G[j],2,Gb[j]);
 148 //printf("Gb\n");
 149 //for(int m=0;m<Gb[j].size();m++) printf("%d ",Gb[j][m]); printf("\n");
 150   }
 151
 152   // alphabet size of each shift register
 153   std::vector<int> bases_x(k);
 154   for(int j=0;j<k ;j++)
 155     bases_x[j] = 1 << max_mem_x[j];
 156 //printf("bases_x\n");
 157 //for(int j=0;j<max_mem_x.size();j++) printf("%d ",max_mem_x[j]); printf("\n");
 158
 159   d_NS.resize(d_I*d_S);
 160   d_OS.resize(d_I*d_S);
 161
 162   std::vector<int> sx(k);
 163   std::vector<int> nsx(k);
 164   std::vector<int> tx(k);
 165   std::vector<std::vector<int> > tb(k);
 166   for(int j=0;j<k;j++)
 167     tb[j].resize(max_mem+1);
 168   std::vector<int> inb(k);
 169   std::vector<int> outb(n);
 170
 171
 172   for(int s=0;s<d_S;s++) {
 173     dec2bases(s,bases_x,sx); // split s into k values, each representing on of the k shift registers
 174 //printf("state = %d \nstates = ",s);
 175 //for(int j=0;j<sx.size();j++) printf("%d ",sx[j]); printf("\n");
 176     for(int i=0;i<d_I;i++) {
 177       dec2base(i,2,inb); // input in binary
 178 //printf("input = %d \ninputs = ",i);
 179 //for(int j=0;j<inb.size();j++) printf("%d ",inb[j]); printf("\n");
 180
 181       // evaluate next state
 182       for(int j=0;j<k;j++)
 183         nsx[j] = (inb[j]*bases_x[j]+sx[j])/2; // next state (for each shift register) MSB first
 184       d_NS[s*d_I+i]=bases2dec(nsx,bases_x); // collect all values into the new state
 185
 186       // evaluate transitions
 187       for(int j=0;j<k;j++)
 188         tx[j] = inb[j]*bases_x[j]+sx[j]; // transition (for each shift register)MSB first
 189       for(int j=0;j<k;j++) {
 190         dec2base(tx[j],2,tb[j]); // transition in binary
 191 //printf("transition = %d \ntransitions = ",tx[j]);
 192 //for(int m=0;m<tb[j].size();m++) printf("%d ",tb[j][m]); printf("\n");
 193       }
 194
 195       // evaluate outputs
 196       for(int nn=0;nn<n;nn++) {
 197         outb[nn] = 0;
 198         for(int j=0;j<k;j++) {
 199           for(int m=0;m<max_mem+1;m++)
 200             outb[nn] = (outb[nn] + Gb[j*n+nn][m]*tb[j][m]) % 2; // careful: polynomial 1+D ir represented as 110, not as 011
 201 //printf("output %d equals %d\n",nn,outb[nn]);
 202         }
 203       }
 204       d_OS[s*d_I+i] = base2dec(outb,2);
 205     }
 206   }
 207
 208   generate_PS_PI();
 209   generate_TM();
 210 }
 211
 212
 213
 214
 215 //######################################################################
 216 //# Automatically generate an FSM specification describing the
 217 //# ISI for a channel
 218 //# of length ch_length and a modulation of size mod_size
 219 //######################################################################
 220 fsm::fsm(int mod_size, int ch_length)
 221 {
 222   d_I=mod_size;
 223   d_S=(int) (pow(1.0*d_I,1.0*ch_length-1)+0.5);
 224   d_O=d_S*d_I;
 225
 226   d_NS.resize(d_I*d_S);
 227   d_OS.resize(d_I*d_S);
 228
 229   for(int s=0;s<d_S;s++) {
 230     for(int i=0;i<d_I;i++) {
 231       int t=i*d_S+s;
 232       d_NS[s*d_I+i] = t/d_I;
 233       d_OS[s*d_I+i] = t;
 234     }
 235   }
 236
 237   generate_PS_PI();
 238   generate_TM();
 239 }
 240
 241
 242 //######################################################################
 243 //# generate the PS and PI tables for later use
 244 //######################################################################
 245 void fsm::generate_PS_PI()
 246 {
 247   d_PS.resize(d_S);
 248   d_PI.resize(d_S);
 249
 250   for(int i=0;i<d_S;i++) {
 251     d_PS[i].resize(d_I*d_S); // max possible size
 252     d_PI[i].resize(d_I*d_S);
 253     int j=0;
 254     for(int ii=0;ii<d_S;ii++) for(int jj=0;jj<d_I;jj++) {
 255       if(d_NS[ii*d_I+jj]!=i) continue;
 256       d_PS[i][j]=ii;
 257       d_PI[i][j]=jj;
 258       j++;
 259     }
 260     d_PS[i].resize(j);
 261     d_PI[i].resize(j);
 262   }
 263 }
 264
 265
 266 //######################################################################
 267 //# generate the termination matrices TMl and TMi for later use
 268 //######################################################################
 269 void fsm::generate_TM()
 270 {
 271   d_TMi.resize(d_S*d_S);
 272   d_TMl.resize(d_S*d_S);
 273
 274   for(int i=0;i<d_S*d_S;i++) {
 275     d_TMi[i] = -1; // no meaning
 276     d_TMl[i] = d_S; //infinity: you need at most S-1 steps
 277     if (i/d_S == i%d_S)
 278       d_TMl[i] = 0;
 279   }
 280
 281   for(int s=0;s<d_S;s++) {
 282     bool done = false;
 283     int attempts = 0;
 284     while (done == false && attempts < d_S-1) {
 285       done = find_es(s);
 286       attempts ++;
 287     }
 288     if (done == false) {
 289       //throw std::runtime_error ("fsm::generate_TM(): FSM appears to be disconnected\n");
 290       printf("fsm::generate_TM(): FSM appears to be disconnected\n");
 291       printf("state %d cannot be reached from all other states\n",s);
 292     }
 293   }
 294 }
 295
 296
 297 // find a path from any state to the ending state "es"
 298 bool fsm::find_es(int es)
 299 {
 300   bool done = true;
 301   for(int s=0;s<d_S;s++) {
 302     if(d_TMl[s*d_S+es] < d_S)
 303       continue;
 304     int minl=d_S;
 305     int mini=-1;
 306     for(int i=0;i<d_I;i++) {
 307       if( 1 + d_TMl[d_NS[s*d_I+i]*d_S+es] < minl) {
 308         minl = 1 + d_TMl[d_NS[s*d_I+i]*d_S+es];
 309         mini = i;
 310       }
 311     }
 312     if (mini != -1) {
 313       d_TMl[s*d_S+es]=minl;
 314       d_TMi[s*d_S+es]=mini;
 315     }
 316     else
 317       done = false;
 318   }
 319   return done;
 320 }
 321
 322
 323
 324
 325
 326 //######################################################################
 327 //#  generate trellis representation of FSM as an SVG file
 328 //######################################################################
 329 void fsm::write_trellis_svg( std::string filename ,int number_stages)
 330 {
 331    std::ofstream trellis_fname (filename.c_str());
 332    if (!trellis_fname) {std::cout << "file not found " << std::endl ; exit(-1);}
 333    const int TRELLIS_Y_OFFSET = 30;
 334    const int TRELLIS_X_OFFSET = 20;
 335    const int STAGE_LABEL_Y_OFFSET = 25;
 336    const int STAGE_LABEL_X_OFFSET = 20;
 337    const int STATE_LABEL_Y_OFFSET = 30;
 338    const int STATE_LABEL_X_OFFSET = 5;
 339    const int STAGE_STATE_OFFSETS = 10;
 340 //   std::cout << "################## BEGIN SVG TRELLIS PIC #####################" << std::endl;
 341    trellis_fname << "<svg viewBox = \"0 0 200 200\" version = \"1.1\">" << std::endl;
 342
 343     for(unsigned int stage_num = 0;stage_num < number_stages;stage_num ++){
 344     // draw states
 345       for (unsigned int state_num = 0;state_num < d_S ; state_num ++ ) {
 346         trellis_fname << "<circle cx = \"" << stage_num * STAGE_STATE_OFFSETS + TRELLIS_X_OFFSET <<
 347         "\" cy = \"" << state_num * STAGE_STATE_OFFSETS + TRELLIS_Y_OFFSET << "\" r = \"1\"/>" << std::endl;
 348       //draw branches
 349         if(stage_num != number_stages-1){
 350           for(unsigned int branch_num = 0;branch_num < d_I; branch_num++){
 351             trellis_fname << "<line x1 =\"" << STAGE_STATE_OFFSETS * stage_num+ TRELLIS_X_OFFSET  << "\" ";
 352             trellis_fname << "y1 =\"" << state_num * STAGE_STATE_OFFSETS + TRELLIS_Y_OFFSET<< "\" ";
 353             trellis_fname << "x2 =\"" <<  STAGE_STATE_OFFSETS *stage_num + STAGE_STATE_OFFSETS+ TRELLIS_X_OFFSET << "\" ";
 354             trellis_fname << "y2 =\"" << d_NS[d_I * state_num + branch_num] * STAGE_STATE_OFFSETS + TRELLIS_Y_OFFSET << "\" ";
 355             trellis_fname << " stroke-dasharray = \"3," <<  branch_num << "\" ";
 356             trellis_fname << " stroke = \"black\" stroke-width = \"0.3\"/>" << std::endl;
 357           }
 358         }
 359       }
 360     }
 361   // label the stages
 362   trellis_fname << "<g font-size = \"4\" font= \"times\" fill = \"black\">" << std::endl;
 363   for(unsigned int stage_num = 0;stage_num < number_stages ;stage_num ++){
 364     trellis_fname << "<text x = \"" << stage_num * STAGE_STATE_OFFSETS + STAGE_LABEL_X_OFFSET <<
 365       "\" y = \""  << STAGE_LABEL_Y_OFFSET  << "\" >" << std::endl;
 366     trellis_fname << stage_num <<  std::endl;
 367     trellis_fname << "</text>" << std::endl;
 368   }
 369   trellis_fname << "</g>" << std::endl;
 370
 371   // label the states
 372   trellis_fname << "<g font-size = \"4\" font= \"times\" fill = \"black\">" << std::endl;
 373   for(unsigned int state_num = 0;state_num < d_S ; state_num ++){
 374     trellis_fname << "<text y = \"" << state_num * STAGE_STATE_OFFSETS + STATE_LABEL_Y_OFFSET <<
 375       "\" x = \""  << STATE_LABEL_X_OFFSET  << "\" >" << std::endl;
 376     trellis_fname << state_num <<  std::endl;
 377     trellis_fname << "</text>" << std::endl;
 378   }
 379   trellis_fname << "</g>" << std::endl;
 380
 381
 382   trellis_fname << "</svg>" << std::endl;
 383 //  std::cout << "################## END SVG TRELLIS PIC ##################### " << std::endl;
 384   trellis_fname.close();
 385 }