3 * Copyright 2004 Free Software Foundation, Inc.
5 * This file is part of GNU Radio
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 3, or (at your option)
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.
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.
27 #include <trellis_siso_combined_f.h>
28 #include <gr_io_signature.h>
33 static const float INF = 1.0e9;
35 trellis_siso_combined_f_sptr
36 trellis_make_siso_combined_f (
43 trellis_siso_type_t SISO_TYPE,
45 const std::vector<float> &TABLE,
46 trellis_metric_type_t TYPE)
48 return trellis_siso_combined_f_sptr (new trellis_siso_combined_f (FSM,K,S0,SK,POSTI,POSTO,SISO_TYPE,D,TABLE,TYPE));
51 trellis_siso_combined_f::trellis_siso_combined_f (
58 trellis_siso_type_t SISO_TYPE,
60 const std::vector<float> &TABLE,
61 trellis_metric_type_t TYPE)
62 : gr_block ("siso_combined_f",
63 gr_make_io_signature (1, -1, sizeof (float)),
64 gr_make_io_signature (1, -1, sizeof (float))),
71 d_SISO_TYPE (SISO_TYPE),
75 //d_alpha(FSM.S()*(K+1)),
76 //d_beta(FSM.S()*(K+1))
79 if (d_POSTI && d_POSTO)
80 multiple = d_FSM.I()+d_FSM.O();
86 throw std::runtime_error ("Not both POSTI and POSTO can be false.");
87 //printf("constructor: Multiple = %d\n",multiple);
88 set_output_multiple (d_K*multiple);
89 //what is the meaning of relative rate for a block with 2 inputs?
90 //set_relative_rate ( multiple / ((double) d_FSM.I()) );
91 // it turns out that the above gives problems in the scheduler, so
92 // let's try (assumption O>I)
93 //set_relative_rate ( multiple / ((double) d_FSM.O()) );
94 // I am tempted to automate like this
96 set_relative_rate ( multiple / ((double) d_D) );
98 set_relative_rate ( multiple / ((double) d_FSM.I()) );
103 trellis_siso_combined_f::forecast (int noutput_items, gr_vector_int &ninput_items_required)
106 if (d_POSTI && d_POSTO)
107 multiple = d_FSM.I()+d_FSM.O();
109 multiple = d_FSM.I();
111 multiple = d_FSM.O();
113 throw std::runtime_error ("Not both POSTI and POSTO can be false.");
114 //printf("forecast: Multiple = %d\n",multiple);
115 assert (noutput_items % (d_K*multiple) == 0);
116 int input_required1 = d_FSM.I() * (noutput_items/multiple) ;
117 int input_required2 = d_D * (noutput_items/multiple) ;
118 //printf("forecast: Output requirements: %d\n",noutput_items);
119 //printf("forecast: Input requirements: %d %d\n",input_required1,input_required2);
120 unsigned ninputs = ninput_items_required.size();
121 assert(ninputs % 2 == 0);
122 for (unsigned int i = 0; i < ninputs/2; i++) {
123 ninput_items_required[2*i] = input_required1;
124 ninput_items_required[2*i+1] = input_required2;
128 inline float min(float a, float b)
130 return a <= b ? a : b;
133 inline float min_star(float a, float b)
135 return (a <= b ? a : b)-log(1+exp(a <= b ? a-b : b-a));
138 void siso_algorithm_combined(int I, int S, int O,
139 const std::vector<int> &NS,
140 const std::vector<int> &OS,
141 const std::vector<int> &PS,
142 const std::vector<int> &PI,
145 bool POSTI, bool POSTO,
146 float (*p2mymin)(float,float),
148 const std::vector<float> &TABLE,
149 trellis_metric_type_t TYPE,
150 const float *priori, const float *observations, float *post//,
151 //std::vector<float> &alpha,
152 //std::vector<float> &beta
156 std::vector<float> alpha(S*(K+1));
157 std::vector<float> beta(S*(K+1));
158 float *prioro = new float[O*K];
161 if(S0<0) { // initial state not specified
162 for(int i=0;i<S;i++) alpha[0*S+i]=0;
165 for(int i=0;i<S;i++) alpha[0*S+i]=INF;
169 for(int k=0;k<K;k++) { // forward recursion
170 calc_metric(O, D, TABLE, &(observations[k*D]), &(prioro[k*O]),TYPE); // calc metrics
172 for(int j=0;j<S;j++) {
174 for(int i=0;i<I;i++) {
176 mm=alpha[k*S+PS[i0]]+priori[k*I+PI[i0]]+prioro[k*O+OS[PS[i0]*I+PI[i0]]];
177 minm=(*p2mymin)(minm,mm);
179 alpha[(k+1)*S+j]=minm;
180 if(minm<norm) norm=minm;
183 alpha[(k+1)*S+j]-=norm; // normalize total metrics so they do not explode
186 if(SK<0) { // final state not specified
187 for(int i=0;i<S;i++) beta[K*S+i]=0;
190 for(int i=0;i<S;i++) beta[K*S+i]=INF;
194 for(int k=K-1;k>=0;k--) { // backward recursion
196 for(int j=0;j<S;j++) {
198 for(int i=0;i<I;i++) {
200 mm=beta[(k+1)*S+NS[i0]]+priori[k*I+i]+prioro[k*O+OS[i0]];
201 minm=(*p2mymin)(minm,mm);
204 if(minm<norm) norm=minm;
207 beta[k*S+j]-=norm; // normalize total metrics so they do not explode
213 for(int k=0;k<K;k++) { // input combining
215 for(int i=0;i<I;i++) {
217 for(int j=0;j<S;j++) {
218 mm=alpha[k*S+j]+prioro[k*O+OS[j*I+i]]+beta[(k+1)*S+NS[j*I+i]];
219 minm=(*p2mymin)(minm,mm);
221 post[k*(I+O)+i]=minm;
222 if(minm<norm) norm=minm;
225 post[k*(I+O)+i]-=norm; // normalize metrics
229 for(int k=0;k<K;k++) { // output combining
231 for(int n=0;n<O;n++) {
233 for(int j=0;j<S;j++) {
234 for(int i=0;i<I;i++) {
235 mm= (n==OS[j*I+i] ? alpha[k*S+j]+priori[k*I+i]+beta[(k+1)*S+NS[j*I+i]] : INF);
236 minm=(*p2mymin)(minm,mm);
239 post[k*(I+O)+I+n]=minm;
240 if(minm<norm) norm=minm;
243 post[k*(I+O)+I+n]-=norm; // normalize metrics
248 for(int k=0;k<K;k++) { // input combining
250 for(int i=0;i<I;i++) {
252 for(int j=0;j<S;j++) {
253 mm=alpha[k*S+j]+prioro[k*O+OS[j*I+i]]+beta[(k+1)*S+NS[j*I+i]];
254 minm=(*p2mymin)(minm,mm);
257 if(minm<norm) norm=minm;
260 post[k*I+i]-=norm; // normalize metrics
265 for(int k=0;k<K;k++) { // output combining
267 for(int n=0;n<O;n++) {
269 for(int j=0;j<S;j++) {
270 for(int i=0;i<I;i++) {
271 mm= (n==OS[j*I+i] ? alpha[k*S+j]+priori[k*I+i]+beta[(k+1)*S+NS[j*I+i]] : INF);
272 minm=(*p2mymin)(minm,mm);
276 if(minm<norm) norm=minm;
279 post[k*O+n]-=norm; // normalize metrics
283 throw std::runtime_error ("Not both POSTI and POSTO can be false.");
295 trellis_siso_combined_f::general_work (int noutput_items,
296 gr_vector_int &ninput_items,
297 gr_vector_const_void_star &input_items,
298 gr_vector_void_star &output_items)
300 assert (input_items.size() == 2*output_items.size());
301 int nstreams = output_items.size();
302 //printf("general_work:Streams: %d\n",nstreams);
304 if (d_POSTI && d_POSTO)
305 multiple = d_FSM.I()+d_FSM.O();
307 multiple = d_FSM.I();
309 multiple = d_FSM.O();
311 throw std::runtime_error ("Not both POSTI and POSTO can be false.");
313 assert (noutput_items % (d_K*multiple) == 0);
314 int nblocks = noutput_items / (d_K*multiple);
315 //printf("general_work:Blocks: %d\n",nblocks);
316 //for(int i=0;i<ninput_items.size();i++)
317 //printf("general_work:Input items available: %d\n",ninput_items[i]);
319 float (*p2min)(float, float) = NULL;
320 if(d_SISO_TYPE == TRELLIS_MIN_SUM)
322 else if(d_SISO_TYPE == TRELLIS_SUM_PRODUCT)
326 for (int m=0;m<nstreams;m++) {
327 const float *in1 = (const float *) input_items[2*m];
328 const float *in2 = (const float *) input_items[2*m+1];
329 float *out = (float *) output_items[m];
330 for (int n=0;n<nblocks;n++) {
331 siso_algorithm_combined(d_FSM.I(),d_FSM.S(),d_FSM.O(),
332 d_FSM.NS(),d_FSM.OS(),d_FSM.PS(),d_FSM.PI(),
337 &(in1[n*d_K*d_FSM.I()]),&(in2[n*d_K*d_D]),
338 &(out[n*d_K*multiple])//,
344 for (unsigned int i = 0; i < input_items.size()/2; i++) {
345 consume(2*i,d_FSM.I() * noutput_items / multiple );
346 consume(2*i+1,d_D * noutput_items / multiple );
349 return noutput_items;