source: trunk/GDE/PHYML20130708/phyml/src/pars.c

Last change on this file was 10307, checked in by aboeckma, 11 years ago

added most recent version of phyml

File size: 33.2 KB
Line 
1/*
2
3PHYML :  a program that  computes maximum likelihood  phyLOGenies from
4DNA or AA homoLOGous sequences
5
6Copyright (C) Stephane Guindon. Oct 2003 onward
7
8All parts of  the source except where indicated  are distributed under
9the GNU public licence.  See http://www.opensource.org for details.
10
11*/
12
13
14#include "pars.h"
15
16
17/*********************************************************/
18
19
20int Pars(t_edge *b, t_tree *tree)
21{
22  int site,n_patterns;
23
24  n_patterns = tree->n_pattern;
25 
26  if(!b)
27    {
28      Post_Order_Pars(tree->a_nodes[0],tree->a_nodes[0]->v[0],tree);
29      if(tree->both_sides == YES) Pre_Order_Pars(tree->a_nodes[0],tree->a_nodes[0]->v[0],tree);
30    }
31
32  if(!b) b = tree->a_nodes[0]->b[0];
33
34  tree->c_pars = 0;
35  For(site,n_patterns)
36    {
37      tree->site_pars[site] = 0;
38      tree->curr_site       = site;
39      tree->site_pars[site] = Pars_Core(b,tree);
40      tree->c_pars += tree->site_pars[site] * tree->data->wght[site];
41      /* printf("\n. site %d pars: %d",site,tree->c_pars); */
42    }
43
44  if(tree->is_mixt_tree) 
45    {
46      MIXT_Pars(b,tree);
47    }
48
49  return tree->c_pars;
50}
51
52/*********************************************************/
53
54void Post_Order_Pars(t_node *a, t_node *d, t_tree *tree)
55{
56  int i,dir;
57
58  dir = -1;
59
60  if(d->tax) return;
61  else
62    {
63      For(i,3)
64        {
65          if(d->v[i] != a)
66            Post_Order_Pars(d,d->v[i],tree);
67          else dir = i;
68        }
69      Get_All_Partial_Pars(tree,d->b[dir],a,d);
70    }
71}
72
73/*********************************************************/
74
75void Pre_Order_Pars(t_node *a, t_node *d, t_tree *tree)
76{
77  int i;
78
79  if(d->tax) return;
80  else
81    {
82      For(i,3)
83        {
84          if(d->v[i] != a)
85            {
86              Get_All_Partial_Pars(tree,d->b[i],d->v[i],d);
87              Pre_Order_Pars(d,d->v[i],tree);
88            }
89        }
90    }
91}
92
93/*********************************************************/
94
95void Get_All_Partial_Pars(t_tree *tree, t_edge *b_fcus, t_node *a, t_node *d)
96{
97  Update_P_Pars(tree,b_fcus,d);
98}
99
100/*********************************************************/
101
102void Site_Pars(t_tree *tree)
103{
104  tree->site_pars[tree->curr_site] = Pars_Core(tree->a_nodes[0]->b[0],tree);
105}
106
107/*********************************************************/
108
109void Init_P_Pars_Tips(t_tree *tree)
110{
111  int curr_site,i,j;
112  short int *state_v;
113  int dim1;
114
115  dim1 = tree->mod->ns;
116
117  state_v = (short int *)mCalloc(tree->mod->ns,sizeof(short int));
118 
119
120
121  For(curr_site,tree->data->crunch_len)
122    {
123      For(i,tree->n_otu)
124        {
125          if(tree->a_nodes[i]->b[0]->rght->tax != 1)
126            {
127              PhyML_Printf("\n== Err in file %s at line %d\n",__FILE__,__LINE__);
128              Exit("\n");           
129            }
130
131          if(tree->io->datatype == NT)
132            {
133              Init_Tips_At_One_Site_Nucleotides_Int(tree->a_nodes[i]->c_seq->state[curr_site],
134                                                    0,
135                                                    state_v);         
136              For(j,tree->mod->ns) tree->a_nodes[i]->b[0]->p_pars_r[curr_site*dim1+j] = MAX_PARS;
137              For(j,tree->mod->ns) if(state_v[j] > 0.5) tree->a_nodes[i]->b[0]->p_pars_r[curr_site*dim1+j] =  0;
138            }
139          else if(tree->io->datatype == AA)
140            {
141              Init_Tips_At_One_Site_AA_Int(tree->a_nodes[i]->c_seq->state[curr_site],
142                                           0,
143                                           state_v);
144              For(j,tree->mod->ns) tree->a_nodes[i]->b[0]->p_pars_r[curr_site*dim1+j] = MAX_PARS;
145              For(j,tree->mod->ns) if(state_v[j] > 0.5) tree->a_nodes[i]->b[0]->p_pars_r[curr_site*dim1+j] =  0;
146            }
147          else if(tree->io->datatype == GENERIC)
148            {
149              Init_Tips_At_One_Site_Generic_Int(tree->a_nodes[i]->c_seq->state+curr_site*tree->mod->io->state_len,
150                                                tree->mod->ns,
151                                                tree->mod->io->state_len,
152                                                0,
153                                                state_v);
154              For(j,tree->mod->ns) tree->a_nodes[i]->b[0]->p_pars_r[curr_site*dim1+j] = MAX_PARS;
155              For(j,tree->mod->ns) if(state_v[j] > 0.5) tree->a_nodes[i]->b[0]->p_pars_r[curr_site*dim1+j] =  0;
156            }
157        }
158    }
159  Free(state_v);
160}
161
162/*********************************************************/
163
164void Init_Ui_Tips(t_tree *tree)
165{ 
166  int curr_site,i,j,br;
167  short int *state_v;
168
169  state_v = (short int *)mCalloc(tree->mod->ns,sizeof(short int));
170 
171  For(curr_site,tree->data->crunch_len)
172    {
173      For(i,tree->n_otu)
174        {
175          if(tree->io->datatype == NT)
176            {
177              if(tree->a_nodes[i]->b[0]->rght->tax != 1)
178                {
179                  PhyML_Printf("\n== Err in file %s at line %d\n",__FILE__,__LINE__);
180                  Exit("\n");       
181                }
182
183              Init_Tips_At_One_Site_Nucleotides_Int(tree->a_nodes[i]->c_seq->state[curr_site],
184                                                    0,
185                                                    state_v);         
186              /* Init_Tips_At_One_Site_Nucleotides_Int(tree->data->c_seq[i]->state[curr_site], */
187              /*                                            0, */
188              /*                                            state_v);          */
189              tree->a_nodes[i]->b[0]->ui_r[curr_site] = 0;
190              For(j,tree->mod->ns) tree->a_nodes[i]->b[0]->ui_r[curr_site] += (unsigned int)(state_v[j] * POW(2,j));
191            }
192          else if(tree->io->datatype == AA)
193            {
194              Init_Tips_At_One_Site_AA_Int(tree->a_nodes[i]->c_seq->state[curr_site],
195                                           0,
196                                           state_v);
197              /* Init_Tips_At_One_Site_AA_Int(tree->data->c_seq[i]->state[curr_site], */
198              /*                                   0, */
199              /*                                   state_v); */
200              tree->a_nodes[i]->b[0]->ui_r[curr_site] = 0;
201              For(j,tree->mod->ns) tree->a_nodes[i]->b[0]->ui_r[curr_site] += (unsigned int)(state_v[j] * POW(2,j));
202            }
203          else if(tree->io->datatype == GENERIC)
204            {
205              Init_Tips_At_One_Site_Generic_Int(tree->a_nodes[i]->c_seq->state+curr_site*tree->mod->io->state_len,
206                                                tree->mod->ns,
207                                                tree->mod->io->state_len,
208                                                0,
209                                                state_v);
210              /* Init_Tips_At_One_Site_Generic_Int(tree->data->c_seq[i]->state+curr_site*tree->mod->io->state_len, */
211              /*                                        tree->mod->ns, */
212              /*                                        tree->mod->io->state_len, */
213              /*                                        0, */
214              /*                                        state_v); */
215              tree->a_nodes[i]->b[0]->ui_r[curr_site] = 0;
216              For(j,tree->mod->ns) tree->a_nodes[i]->b[0]->ui_r[curr_site] += (unsigned int)(state_v[j] * POW(2,j));
217            }
218        }
219    }
220
221 
222  For(br,2*tree->n_otu-3)
223    {
224      For(curr_site,tree->data->crunch_len)
225        {
226          tree->a_edges[br]->pars_r[curr_site] = 0;
227          tree->a_edges[br]->pars_l[curr_site] = 0;
228        }
229    }
230
231
232  Free(state_v);
233}
234
235/*********************************************************/
236
237void Update_P_Pars(t_tree *tree, t_edge *b_fcus, t_node *n)
238{
239/* 
240           |
241           |<- b_cus
242           |
243           n
244          / \
245         /   \
246        /     \
247*/
248
249  int i,j;
250  int site;
251  unsigned int *ui, *ui_v1, *ui_v2;
252  int *p_pars_v1, *p_pars_v2, *p_pars;
253  int *pars, *pars_v1, *pars_v2;
254  int n_patterns;
255  int min_v1,min_v2;
256  int v;
257  int dim1;
258
259  if((tree->io->do_alias_subpatt == YES) && 
260     (tree->update_alias_subpatt == YES))
261    Alias_One_Subpatt((n==b_fcus->left)?(b_fcus->rght):(b_fcus->left),n,tree);
262
263  if(n->tax) return;
264
265  dim1 = tree->mod->ns;
266  ui = ui_v1 = ui_v2 = NULL;
267  p_pars = p_pars_v1 = p_pars_v2 = NULL;
268  pars = pars_v1 = pars_v2 = NULL;
269
270  n_patterns = tree->n_pattern;
271
272
273  if(n == b_fcus->left)
274    {       
275      ui = b_fcus->ui_l;
276
277      pars = b_fcus->pars_l;
278      p_pars = b_fcus->p_pars_l;
279
280      ui_v1 = 
281      (n == n->b[b_fcus->l_v1]->left)?
282      (n->b[b_fcus->l_v1]->ui_r):
283      (n->b[b_fcus->l_v1]->ui_l);
284
285      ui_v2 = 
286      (n == n->b[b_fcus->l_v2]->left)?
287      (n->b[b_fcus->l_v2]->ui_r):
288      (n->b[b_fcus->l_v2]->ui_l);
289
290      p_pars_v1 = 
291      (n == n->b[b_fcus->l_v1]->left)?
292      (n->b[b_fcus->l_v1]->p_pars_r):
293      (n->b[b_fcus->l_v1]->p_pars_l);
294
295      p_pars_v2 = 
296      (n == n->b[b_fcus->l_v2]->left)?
297      (n->b[b_fcus->l_v2]->p_pars_r):
298      (n->b[b_fcus->l_v2]->p_pars_l);
299
300      pars_v1 = 
301      (n == n->b[b_fcus->l_v1]->left)?
302      (n->b[b_fcus->l_v1]->pars_r):
303      (n->b[b_fcus->l_v1]->pars_l);
304
305      pars_v2 = 
306      (n == n->b[b_fcus->l_v2]->left)?
307      (n->b[b_fcus->l_v2]->pars_r):
308      (n->b[b_fcus->l_v2]->pars_l);
309    }
310  else
311    {
312      ui = b_fcus->ui_r;
313     
314      pars = b_fcus->pars_r;
315      p_pars = b_fcus->p_pars_r;
316
317      ui_v1 = 
318      (n == n->b[b_fcus->r_v1]->left)?
319      (n->b[b_fcus->r_v1]->ui_r):
320      (n->b[b_fcus->r_v1]->ui_l);
321
322      ui_v2 = 
323      (n == n->b[b_fcus->r_v2]->left)?
324      (n->b[b_fcus->r_v2]->ui_r):
325      (n->b[b_fcus->r_v2]->ui_l);
326
327      p_pars_v1 = 
328      (n == n->b[b_fcus->r_v1]->left)?
329      (n->b[b_fcus->r_v1]->p_pars_r):
330      (n->b[b_fcus->r_v1]->p_pars_l);
331
332      p_pars_v2 = 
333      (n == n->b[b_fcus->r_v2]->left)?
334      (n->b[b_fcus->r_v2]->p_pars_r):
335      (n->b[b_fcus->r_v2]->p_pars_l);
336
337      pars_v1 = 
338      (n == n->b[b_fcus->r_v1]->left)?
339      (n->b[b_fcus->r_v1]->pars_r):
340      (n->b[b_fcus->r_v1]->pars_l);
341
342      pars_v2 = 
343      (n == n->b[b_fcus->r_v2]->left)?
344      (n->b[b_fcus->r_v2]->pars_r):
345      (n->b[b_fcus->r_v2]->pars_l);
346    }
347
348
349  if(tree->mod->s_opt->general_pars)
350    {
351      For(site,n_patterns)
352        {
353          For(i,tree->mod->ns)
354            {
355              min_v1 = MAX_PARS;
356              For(j,tree->mod->ns)
357                {
358                  v = p_pars_v1[site*dim1+j] + tree->step_mat[i*tree->mod->ns+j]; 
359                  if(v < min_v1) min_v1 = v;
360                }
361             
362              min_v2 = MAX_PARS;
363              For(j,tree->mod->ns)
364                {
365                  v = p_pars_v2[site*dim1+j] + tree->step_mat[i*tree->mod->ns+j]; 
366                  if(v < min_v2) min_v2 = v;
367                }             
368              p_pars[site*dim1+i] = min_v1 + min_v2;
369            }     
370        }
371    }
372  else
373    {
374      For(site,n_patterns)
375        {
376          pars[site] = pars_v1[site] + pars_v2[site];
377
378          ui[site] = ui_v1[site] & ui_v2[site];
379
380          if(!ui[site])
381            {
382              pars[site]++;
383              ui[site] = ui_v1[site] | ui_v2[site];
384            }
385
386        }
387    }
388}
389
390/*********************************************************/
391
392int Pars_Core(t_edge *b, t_tree *tree)
393{
394  int site;
395  int i,j;
396  int site_pars;
397  int min_l,min_r;
398  int v;
399  int dim1;
400
401  dim1 = tree->mod->ns;
402  site = tree->curr_site;
403  site_pars = MAX_PARS;
404
405  if(tree->mod->s_opt->general_pars)
406    {
407      For(i,tree->mod->ns)
408        {
409          min_l = MAX_PARS;
410          For(j,tree->mod->ns)
411            {
412              v = b->p_pars_l[site*dim1+j] + tree->step_mat[i*tree->mod->ns+j]; 
413              if(v < min_l) min_l = v;
414            }
415
416          min_r = MAX_PARS;
417          For(j,tree->mod->ns)
418            {
419              v = b->p_pars_r[site*dim1+j] + tree->step_mat[i*tree->mod->ns+j]; 
420              if(v < min_r) min_r = v;
421            }
422         
423          if((min_l + min_r) < site_pars) site_pars = min_l + min_r;
424        }
425    }
426  else
427    {
428      site_pars = b->pars_l[site] + b->pars_r[site];     
429      if(!(b->ui_l[site] & b->ui_r[site])) site_pars++;
430    }
431
432  return site_pars;
433}
434
435/*********************************************************/
436/* Is there one or more parsimoniy step(s) along this t_edge ?
437   0 -> NO; 1 -> YES
438*/
439int One_Pars_Step(t_edge *b,t_tree *tree)
440{
441  int site;
442  int init_general_pars;
443
444  init_general_pars = tree->mod->s_opt->general_pars;
445 
446  tree->mod->s_opt->general_pars = 0;
447  Set_Both_Sides(YES,tree);
448  Pars(NULL,tree);
449
450  For(site,tree->n_pattern)
451    {
452      if(!(b->ui_l[site] & b->ui_r[site])) break;
453    }
454  tree->mod->s_opt->general_pars = init_general_pars;
455  if(site == tree->n_pattern) return 0;
456  else                       
457    { 
458      PhyML_Printf("\n. One parsimony step ocurred at site %4d",site);
459      return 1;
460    }
461}
462
463/*********************************************************/
464int Pars_At_Given_Edge(t_edge *b, t_tree *tree)
465{
466  int site,n_patterns;
467
468/*   n_patterns = (int)FLOOR(tree->n_pattern*tree->prop_of_sites_to_consider); */
469  n_patterns = tree->n_pattern;
470
471  tree->c_pars = .0;
472  For(site,n_patterns)
473    {
474      tree->site_pars[site] = 0;
475      tree->curr_site = site;
476      tree->site_pars[site] = Pars_Core(b,tree);
477      tree->c_pars += tree->site_pars[site] * tree->data->wght[site];
478    }
479  return tree->c_pars;
480}
481
482/*********************************************************/
483
484int Update_Pars_At_Given_Edge(t_edge *b_fcus, t_tree *tree)
485{
486  Update_P_Pars(tree,b_fcus,b_fcus->left);
487  Update_P_Pars(tree,b_fcus,b_fcus->rght);
488  tree->c_pars = Pars(b_fcus,tree);
489  return tree->c_pars;
490}
491
492/*********************************************************/
493
494void Get_Step_Mat(t_tree *tree)
495{
496  int i;
497
498  if(tree->io->datatype == AA)
499    {
500      tree->step_mat[ 0*tree->mod->ns+ 0] =    0 ;
501      tree->step_mat[ 0*tree->mod->ns+ 1] =    3 ;
502      tree->step_mat[ 0*tree->mod->ns+ 2] =    3 ;
503      tree->step_mat[ 0*tree->mod->ns+ 3] =    2 ;
504      tree->step_mat[ 0*tree->mod->ns+ 4] =    3 ;
505      tree->step_mat[ 0*tree->mod->ns+ 5] =    3 ;
506      tree->step_mat[ 0*tree->mod->ns+ 6] =    2 ;
507      tree->step_mat[ 0*tree->mod->ns+ 7] =    2 ;
508      tree->step_mat[ 0*tree->mod->ns+ 8] =    3 ;
509      tree->step_mat[ 0*tree->mod->ns+ 9] =    3 ;
510      tree->step_mat[ 0*tree->mod->ns+10] =    3 ;
511      tree->step_mat[ 0*tree->mod->ns+11] =    3 ;
512      tree->step_mat[ 0*tree->mod->ns+12] =    3 ;
513      tree->step_mat[ 0*tree->mod->ns+13] =    3 ;
514      tree->step_mat[ 0*tree->mod->ns+14] =    2 ;
515      tree->step_mat[ 0*tree->mod->ns+15] =    2 ;
516      tree->step_mat[ 0*tree->mod->ns+16] =    2 ;
517      tree->step_mat[ 0*tree->mod->ns+17] =    3 ;
518      tree->step_mat[ 0*tree->mod->ns+18] =    3 ;
519      tree->step_mat[ 0*tree->mod->ns+19] =    2 ;
520      tree->step_mat[ 1*tree->mod->ns+ 0] =    3 ;
521      tree->step_mat[ 1*tree->mod->ns+ 1] =    0 ;
522      tree->step_mat[ 1*tree->mod->ns+ 2] =    2 ;
523      tree->step_mat[ 1*tree->mod->ns+ 3] =    3 ;
524      tree->step_mat[ 1*tree->mod->ns+ 4] =    2 ;
525      tree->step_mat[ 1*tree->mod->ns+ 5] =    2 ;
526      tree->step_mat[ 1*tree->mod->ns+ 6] =    3 ;
527      tree->step_mat[ 1*tree->mod->ns+ 7] =    2 ;
528      tree->step_mat[ 1*tree->mod->ns+ 8] =    2 ;
529      tree->step_mat[ 1*tree->mod->ns+ 9] =    2 ;
530      tree->step_mat[ 1*tree->mod->ns+10] =    2 ;
531      tree->step_mat[ 1*tree->mod->ns+11] =    2 ;
532      tree->step_mat[ 1*tree->mod->ns+12] =    2 ;
533      tree->step_mat[ 1*tree->mod->ns+13] =    3 ;
534      tree->step_mat[ 1*tree->mod->ns+14] =    2 ;
535      tree->step_mat[ 1*tree->mod->ns+15] =    2 ;
536      tree->step_mat[ 1*tree->mod->ns+16] =    2 ;
537      tree->step_mat[ 1*tree->mod->ns+17] =    2 ;
538      tree->step_mat[ 1*tree->mod->ns+18] =    3 ;
539      tree->step_mat[ 1*tree->mod->ns+19] =    3 ;
540      tree->step_mat[ 2*tree->mod->ns+ 0] =    3 ;
541      tree->step_mat[ 2*tree->mod->ns+ 1] =    2 ;
542      tree->step_mat[ 2*tree->mod->ns+ 2] =    0 ;
543      tree->step_mat[ 2*tree->mod->ns+ 3] =    2 ;
544      tree->step_mat[ 2*tree->mod->ns+ 4] =    2 ;
545      tree->step_mat[ 2*tree->mod->ns+ 5] =    2 ;
546      tree->step_mat[ 2*tree->mod->ns+ 6] =    2 ;
547      tree->step_mat[ 2*tree->mod->ns+ 7] =    3 ;
548      tree->step_mat[ 2*tree->mod->ns+ 8] =    2 ;
549      tree->step_mat[ 2*tree->mod->ns+ 9] =    2 ;
550      tree->step_mat[ 2*tree->mod->ns+10] =    3 ;
551      tree->step_mat[ 2*tree->mod->ns+11] =    1 ;
552      tree->step_mat[ 2*tree->mod->ns+12] =    2 ;
553      tree->step_mat[ 2*tree->mod->ns+13] =    2 ;
554      tree->step_mat[ 2*tree->mod->ns+14] =    3 ;
555      tree->step_mat[ 2*tree->mod->ns+15] =    2 ;
556      tree->step_mat[ 2*tree->mod->ns+16] =    2 ;
557      tree->step_mat[ 2*tree->mod->ns+17] =    3 ;
558      tree->step_mat[ 2*tree->mod->ns+18] =    2 ;
559      tree->step_mat[ 2*tree->mod->ns+19] =    3 ;
560      tree->step_mat[ 3*tree->mod->ns+ 0] =    2 ;
561      tree->step_mat[ 3*tree->mod->ns+ 1] =    3 ;
562      tree->step_mat[ 3*tree->mod->ns+ 2] =    2 ;
563      tree->step_mat[ 3*tree->mod->ns+ 3] =    0 ;
564      tree->step_mat[ 3*tree->mod->ns+ 4] =    2 ;
565      tree->step_mat[ 3*tree->mod->ns+ 5] =    2 ;
566      tree->step_mat[ 3*tree->mod->ns+ 6] =    1 ;
567      tree->step_mat[ 3*tree->mod->ns+ 7] =    2 ;
568      tree->step_mat[ 3*tree->mod->ns+ 8] =    2 ;
569      tree->step_mat[ 3*tree->mod->ns+ 9] =    3 ;
570      tree->step_mat[ 3*tree->mod->ns+10] =    3 ;
571      tree->step_mat[ 3*tree->mod->ns+11] =    2 ;
572      tree->step_mat[ 3*tree->mod->ns+12] =    3 ;
573      tree->step_mat[ 3*tree->mod->ns+13] =    2 ;
574      tree->step_mat[ 3*tree->mod->ns+14] =    3 ;
575      tree->step_mat[ 3*tree->mod->ns+15] =    3 ;
576      tree->step_mat[ 3*tree->mod->ns+16] =    3 ;
577      tree->step_mat[ 3*tree->mod->ns+17] =    3 ;
578      tree->step_mat[ 3*tree->mod->ns+18] =    2 ;
579      tree->step_mat[ 3*tree->mod->ns+19] =    2 ;
580      tree->step_mat[ 4*tree->mod->ns+ 0] =    3 ;
581      tree->step_mat[ 4*tree->mod->ns+ 1] =    2 ;
582      tree->step_mat[ 4*tree->mod->ns+ 2] =    2 ;
583      tree->step_mat[ 4*tree->mod->ns+ 3] =    2 ;
584      tree->step_mat[ 4*tree->mod->ns+ 4] =    0 ;
585      tree->step_mat[ 4*tree->mod->ns+ 5] =    3 ;
586      tree->step_mat[ 4*tree->mod->ns+ 6] =    3 ;
587      tree->step_mat[ 4*tree->mod->ns+ 7] =    2 ;
588      tree->step_mat[ 4*tree->mod->ns+ 8] =    2 ;
589      tree->step_mat[ 4*tree->mod->ns+ 9] =    3 ;
590      tree->step_mat[ 4*tree->mod->ns+10] =    2 ;
591      tree->step_mat[ 4*tree->mod->ns+11] =    3 ;
592      tree->step_mat[ 4*tree->mod->ns+12] =    3 ;
593      tree->step_mat[ 4*tree->mod->ns+13] =    2 ;
594      tree->step_mat[ 4*tree->mod->ns+14] =    3 ;
595      tree->step_mat[ 4*tree->mod->ns+15] =    2 ;
596      tree->step_mat[ 4*tree->mod->ns+16] =    3 ;
597      tree->step_mat[ 4*tree->mod->ns+17] =    1 ;
598      tree->step_mat[ 4*tree->mod->ns+18] =    2 ;
599      tree->step_mat[ 4*tree->mod->ns+19] =    3 ;
600      tree->step_mat[ 5*tree->mod->ns+ 0] =    3 ;
601      tree->step_mat[ 5*tree->mod->ns+ 1] =    2 ;
602      tree->step_mat[ 5*tree->mod->ns+ 2] =    2 ;
603      tree->step_mat[ 5*tree->mod->ns+ 3] =    2 ;
604      tree->step_mat[ 5*tree->mod->ns+ 4] =    3 ;
605      tree->step_mat[ 5*tree->mod->ns+ 5] =    0 ;
606      tree->step_mat[ 5*tree->mod->ns+ 6] =    2 ;
607      tree->step_mat[ 5*tree->mod->ns+ 7] =    3 ;
608      tree->step_mat[ 5*tree->mod->ns+ 8] =    1 ;
609      tree->step_mat[ 5*tree->mod->ns+ 9] =    3 ;
610      tree->step_mat[ 5*tree->mod->ns+10] =    2 ;
611      tree->step_mat[ 5*tree->mod->ns+11] =    2 ;
612      tree->step_mat[ 5*tree->mod->ns+12] =    2 ;
613      tree->step_mat[ 5*tree->mod->ns+13] =    3 ;
614      tree->step_mat[ 5*tree->mod->ns+14] =    2 ;
615      tree->step_mat[ 5*tree->mod->ns+15] =    3 ;
616      tree->step_mat[ 5*tree->mod->ns+16] =    3 ;
617      tree->step_mat[ 5*tree->mod->ns+17] =    2 ;
618      tree->step_mat[ 5*tree->mod->ns+18] =    2 ;
619      tree->step_mat[ 5*tree->mod->ns+19] =    3 ;
620      tree->step_mat[ 6*tree->mod->ns+ 0] =    2 ;
621      tree->step_mat[ 6*tree->mod->ns+ 1] =    3 ;
622      tree->step_mat[ 6*tree->mod->ns+ 2] =    2 ;
623      tree->step_mat[ 6*tree->mod->ns+ 3] =    1 ;
624      tree->step_mat[ 6*tree->mod->ns+ 4] =    3 ;
625      tree->step_mat[ 6*tree->mod->ns+ 5] =    2 ;
626      tree->step_mat[ 6*tree->mod->ns+ 6] =    0 ;
627      tree->step_mat[ 6*tree->mod->ns+ 7] =    2 ;
628      tree->step_mat[ 6*tree->mod->ns+ 8] =    2 ;
629      tree->step_mat[ 6*tree->mod->ns+ 9] =    3 ;
630      tree->step_mat[ 6*tree->mod->ns+10] =    3 ;
631      tree->step_mat[ 6*tree->mod->ns+11] =    2 ;
632      tree->step_mat[ 6*tree->mod->ns+12] =    2 ;
633      tree->step_mat[ 6*tree->mod->ns+13] =    3 ;
634      tree->step_mat[ 6*tree->mod->ns+14] =    3 ;
635      tree->step_mat[ 6*tree->mod->ns+15] =    3 ;
636      tree->step_mat[ 6*tree->mod->ns+16] =    3 ;
637      tree->step_mat[ 6*tree->mod->ns+17] =    2 ;
638      tree->step_mat[ 6*tree->mod->ns+18] =    2 ;
639      tree->step_mat[ 6*tree->mod->ns+19] =    2 ;
640      tree->step_mat[ 7*tree->mod->ns+ 0] =    2 ;
641      tree->step_mat[ 7*tree->mod->ns+ 1] =    2 ;
642      tree->step_mat[ 7*tree->mod->ns+ 2] =    3 ;
643      tree->step_mat[ 7*tree->mod->ns+ 3] =    2 ;
644      tree->step_mat[ 7*tree->mod->ns+ 4] =    2 ;
645      tree->step_mat[ 7*tree->mod->ns+ 5] =    3 ;
646      tree->step_mat[ 7*tree->mod->ns+ 6] =    2 ;
647      tree->step_mat[ 7*tree->mod->ns+ 7] =    0 ;
648      tree->step_mat[ 7*tree->mod->ns+ 8] =    3 ;
649      tree->step_mat[ 7*tree->mod->ns+ 9] =    3 ;
650      tree->step_mat[ 7*tree->mod->ns+10] =    3 ;
651      tree->step_mat[ 7*tree->mod->ns+11] =    3 ;
652      tree->step_mat[ 7*tree->mod->ns+12] =    3 ;
653      tree->step_mat[ 7*tree->mod->ns+13] =    3 ;
654      tree->step_mat[ 7*tree->mod->ns+14] =    3 ;
655      tree->step_mat[ 7*tree->mod->ns+15] =    2 ;
656      tree->step_mat[ 7*tree->mod->ns+16] =    3 ;
657      tree->step_mat[ 7*tree->mod->ns+17] =    2 ;
658      tree->step_mat[ 7*tree->mod->ns+18] =    3 ;
659      tree->step_mat[ 7*tree->mod->ns+19] =    2 ;
660      tree->step_mat[ 8*tree->mod->ns+ 0] =    3 ;
661      tree->step_mat[ 8*tree->mod->ns+ 1] =    2 ;
662      tree->step_mat[ 8*tree->mod->ns+ 2] =    2 ;
663      tree->step_mat[ 8*tree->mod->ns+ 3] =    2 ;
664      tree->step_mat[ 8*tree->mod->ns+ 4] =    2 ;
665      tree->step_mat[ 8*tree->mod->ns+ 5] =    1 ;
666      tree->step_mat[ 8*tree->mod->ns+ 6] =    2 ;
667      tree->step_mat[ 8*tree->mod->ns+ 7] =    3 ;
668      tree->step_mat[ 8*tree->mod->ns+ 8] =    0 ;
669      tree->step_mat[ 8*tree->mod->ns+ 9] =    3 ;
670      tree->step_mat[ 8*tree->mod->ns+10] =    2 ;
671      tree->step_mat[ 8*tree->mod->ns+11] =    2 ;
672      tree->step_mat[ 8*tree->mod->ns+12] =    3 ;
673      tree->step_mat[ 8*tree->mod->ns+13] =    2 ;
674      tree->step_mat[ 8*tree->mod->ns+14] =    2 ;
675      tree->step_mat[ 8*tree->mod->ns+15] =    3 ;
676      tree->step_mat[ 8*tree->mod->ns+16] =    3 ;
677      tree->step_mat[ 8*tree->mod->ns+17] =    3 ;
678      tree->step_mat[ 8*tree->mod->ns+18] =    2 ;
679      tree->step_mat[ 8*tree->mod->ns+19] =    3 ;
680      tree->step_mat[ 9*tree->mod->ns+ 0] =    3 ;
681      tree->step_mat[ 9*tree->mod->ns+ 1] =    2 ;
682      tree->step_mat[ 9*tree->mod->ns+ 2] =    2 ;
683      tree->step_mat[ 9*tree->mod->ns+ 3] =    3 ;
684      tree->step_mat[ 9*tree->mod->ns+ 4] =    3 ;
685      tree->step_mat[ 9*tree->mod->ns+ 5] =    3 ;
686      tree->step_mat[ 9*tree->mod->ns+ 6] =    3 ;
687      tree->step_mat[ 9*tree->mod->ns+ 7] =    3 ;
688      tree->step_mat[ 9*tree->mod->ns+ 8] =    3 ;
689      tree->step_mat[ 9*tree->mod->ns+ 9] =    0 ;
690      tree->step_mat[ 9*tree->mod->ns+10] =    2 ;
691      tree->step_mat[ 9*tree->mod->ns+11] =    2 ;
692      tree->step_mat[ 9*tree->mod->ns+12] =    1 ;
693      tree->step_mat[ 9*tree->mod->ns+13] =    2 ;
694      tree->step_mat[ 9*tree->mod->ns+14] =    3 ;
695      tree->step_mat[ 9*tree->mod->ns+15] =    2 ;
696      tree->step_mat[ 9*tree->mod->ns+16] =    2 ;
697      tree->step_mat[ 9*tree->mod->ns+17] =    3 ;
698      tree->step_mat[ 9*tree->mod->ns+18] =    3 ;
699      tree->step_mat[ 9*tree->mod->ns+19] =    2 ;
700      tree->step_mat[10*tree->mod->ns+ 0] =    3 ;
701      tree->step_mat[10*tree->mod->ns+ 1] =    2 ;
702      tree->step_mat[10*tree->mod->ns+ 2] =    3 ;
703      tree->step_mat[10*tree->mod->ns+ 3] =    3 ;
704      tree->step_mat[10*tree->mod->ns+ 4] =    2 ;
705      tree->step_mat[10*tree->mod->ns+ 5] =    2 ;
706      tree->step_mat[10*tree->mod->ns+ 6] =    3 ;
707      tree->step_mat[10*tree->mod->ns+ 7] =    3 ;
708      tree->step_mat[10*tree->mod->ns+ 8] =    2 ;
709      tree->step_mat[10*tree->mod->ns+ 9] =    2 ;
710      tree->step_mat[10*tree->mod->ns+10] =    0 ;
711      tree->step_mat[10*tree->mod->ns+11] =    3 ;
712      tree->step_mat[10*tree->mod->ns+12] =    2 ;
713      tree->step_mat[10*tree->mod->ns+13] =    2 ;
714      tree->step_mat[10*tree->mod->ns+14] =    2 ;
715      tree->step_mat[10*tree->mod->ns+15] =    3 ;
716      tree->step_mat[10*tree->mod->ns+16] =    3 ;
717      tree->step_mat[10*tree->mod->ns+17] =    2 ;
718      tree->step_mat[10*tree->mod->ns+18] =    2 ;
719      tree->step_mat[10*tree->mod->ns+19] =    2 ;
720      tree->step_mat[11*tree->mod->ns+ 0] =    3 ;
721      tree->step_mat[11*tree->mod->ns+ 1] =    2 ;
722      tree->step_mat[11*tree->mod->ns+ 2] =    1 ;
723      tree->step_mat[11*tree->mod->ns+ 3] =    2 ;
724      tree->step_mat[11*tree->mod->ns+ 4] =    3 ;
725      tree->step_mat[11*tree->mod->ns+ 5] =    2 ;
726      tree->step_mat[11*tree->mod->ns+ 6] =    2 ;
727      tree->step_mat[11*tree->mod->ns+ 7] =    3 ;
728      tree->step_mat[11*tree->mod->ns+ 8] =    2 ;
729      tree->step_mat[11*tree->mod->ns+ 9] =    2 ;
730      tree->step_mat[11*tree->mod->ns+10] =    3 ;
731      tree->step_mat[11*tree->mod->ns+11] =    0 ;
732      tree->step_mat[11*tree->mod->ns+12] =    2 ;
733      tree->step_mat[11*tree->mod->ns+13] =    3 ;
734      tree->step_mat[11*tree->mod->ns+14] =    3 ;
735      tree->step_mat[11*tree->mod->ns+15] =    2 ;
736      tree->step_mat[11*tree->mod->ns+16] =    2 ;
737      tree->step_mat[11*tree->mod->ns+17] =    2 ;
738      tree->step_mat[11*tree->mod->ns+18] =    2 ;
739      tree->step_mat[11*tree->mod->ns+19] =    3 ;
740      tree->step_mat[12*tree->mod->ns+ 0] =    3 ;
741      tree->step_mat[12*tree->mod->ns+ 1] =    2 ;
742      tree->step_mat[12*tree->mod->ns+ 2] =    2 ;
743      tree->step_mat[12*tree->mod->ns+ 3] =    3 ;
744      tree->step_mat[12*tree->mod->ns+ 4] =    3 ;
745      tree->step_mat[12*tree->mod->ns+ 5] =    2 ;
746      tree->step_mat[12*tree->mod->ns+ 6] =    2 ;
747      tree->step_mat[12*tree->mod->ns+ 7] =    3 ;
748      tree->step_mat[12*tree->mod->ns+ 8] =    3 ;
749      tree->step_mat[12*tree->mod->ns+ 9] =    1 ;
750      tree->step_mat[12*tree->mod->ns+10] =    2 ;
751      tree->step_mat[12*tree->mod->ns+11] =    2 ;
752      tree->step_mat[12*tree->mod->ns+12] =    0 ;
753      tree->step_mat[12*tree->mod->ns+13] =    2 ;
754      tree->step_mat[12*tree->mod->ns+14] =    3 ;
755      tree->step_mat[12*tree->mod->ns+15] =    2 ;
756      tree->step_mat[12*tree->mod->ns+16] =    2 ;
757      tree->step_mat[12*tree->mod->ns+17] =    2 ;
758      tree->step_mat[12*tree->mod->ns+18] =    3 ;
759      tree->step_mat[12*tree->mod->ns+19] =    2 ;
760      tree->step_mat[13*tree->mod->ns+ 0] =    3 ;
761      tree->step_mat[13*tree->mod->ns+ 1] =    3 ;
762      tree->step_mat[13*tree->mod->ns+ 2] =    2 ;
763      tree->step_mat[13*tree->mod->ns+ 3] =    2 ;
764      tree->step_mat[13*tree->mod->ns+ 4] =    2 ;
765      tree->step_mat[13*tree->mod->ns+ 5] =    3 ;
766      tree->step_mat[13*tree->mod->ns+ 6] =    3 ;
767      tree->step_mat[13*tree->mod->ns+ 7] =    3 ;
768      tree->step_mat[13*tree->mod->ns+ 8] =    2 ;
769      tree->step_mat[13*tree->mod->ns+ 9] =    2 ;
770      tree->step_mat[13*tree->mod->ns+10] =    2 ;
771      tree->step_mat[13*tree->mod->ns+11] =    3 ;
772      tree->step_mat[13*tree->mod->ns+12] =    2 ;
773      tree->step_mat[13*tree->mod->ns+13] =    0 ;
774      tree->step_mat[13*tree->mod->ns+14] =    3 ;
775      tree->step_mat[13*tree->mod->ns+15] =    2 ;
776      tree->step_mat[13*tree->mod->ns+16] =    3 ;
777      tree->step_mat[13*tree->mod->ns+17] =    2 ;
778      tree->step_mat[13*tree->mod->ns+18] =    2 ;
779      tree->step_mat[13*tree->mod->ns+19] =    2 ;
780      tree->step_mat[14*tree->mod->ns+ 0] =    2 ;
781      tree->step_mat[14*tree->mod->ns+ 1] =    2 ;
782      tree->step_mat[14*tree->mod->ns+ 2] =    3 ;
783      tree->step_mat[14*tree->mod->ns+ 3] =    3 ;
784      tree->step_mat[14*tree->mod->ns+ 4] =    3 ;
785      tree->step_mat[14*tree->mod->ns+ 5] =    2 ;
786      tree->step_mat[14*tree->mod->ns+ 6] =    3 ;
787      tree->step_mat[14*tree->mod->ns+ 7] =    3 ;
788      tree->step_mat[14*tree->mod->ns+ 8] =    2 ;
789      tree->step_mat[14*tree->mod->ns+ 9] =    3 ;
790      tree->step_mat[14*tree->mod->ns+10] =    2 ;
791      tree->step_mat[14*tree->mod->ns+11] =    3 ;
792      tree->step_mat[14*tree->mod->ns+12] =    3 ;
793      tree->step_mat[14*tree->mod->ns+13] =    3 ;
794      tree->step_mat[14*tree->mod->ns+14] =    0 ;
795      tree->step_mat[14*tree->mod->ns+15] =    2 ;
796      tree->step_mat[14*tree->mod->ns+16] =    2 ;
797      tree->step_mat[14*tree->mod->ns+17] =    3 ;
798      tree->step_mat[14*tree->mod->ns+18] =    3 ;
799      tree->step_mat[14*tree->mod->ns+19] =    3 ;
800      tree->step_mat[15*tree->mod->ns+ 0] =    2 ;
801      tree->step_mat[15*tree->mod->ns+ 1] =    2 ;
802      tree->step_mat[15*tree->mod->ns+ 2] =    2 ;
803      tree->step_mat[15*tree->mod->ns+ 3] =    3 ;
804      tree->step_mat[15*tree->mod->ns+ 4] =    2 ;
805      tree->step_mat[15*tree->mod->ns+ 5] =    3 ;
806      tree->step_mat[15*tree->mod->ns+ 6] =    3 ;
807      tree->step_mat[15*tree->mod->ns+ 7] =    2 ;
808      tree->step_mat[15*tree->mod->ns+ 8] =    3 ;
809      tree->step_mat[15*tree->mod->ns+ 9] =    2 ;
810      tree->step_mat[15*tree->mod->ns+10] =    3 ;
811      tree->step_mat[15*tree->mod->ns+11] =    2 ;
812      tree->step_mat[15*tree->mod->ns+12] =    2 ;
813      tree->step_mat[15*tree->mod->ns+13] =    2 ;
814      tree->step_mat[15*tree->mod->ns+14] =    2 ;
815      tree->step_mat[15*tree->mod->ns+15] =    0 ;
816      tree->step_mat[15*tree->mod->ns+16] =    2 ;
817      tree->step_mat[15*tree->mod->ns+17] =    2 ;
818      tree->step_mat[15*tree->mod->ns+18] =    2 ;
819      tree->step_mat[15*tree->mod->ns+19] =    3 ;
820      tree->step_mat[16*tree->mod->ns+ 0] =    2 ;
821      tree->step_mat[16*tree->mod->ns+ 1] =    2 ;
822      tree->step_mat[16*tree->mod->ns+ 2] =    2 ;
823      tree->step_mat[16*tree->mod->ns+ 3] =    3 ;
824      tree->step_mat[16*tree->mod->ns+ 4] =    3 ;
825      tree->step_mat[16*tree->mod->ns+ 5] =    3 ;
826      tree->step_mat[16*tree->mod->ns+ 6] =    3 ;
827      tree->step_mat[16*tree->mod->ns+ 7] =    3 ;
828      tree->step_mat[16*tree->mod->ns+ 8] =    3 ;
829      tree->step_mat[16*tree->mod->ns+ 9] =    2 ;
830      tree->step_mat[16*tree->mod->ns+10] =    3 ;
831      tree->step_mat[16*tree->mod->ns+11] =    2 ;
832      tree->step_mat[16*tree->mod->ns+12] =    2 ;
833      tree->step_mat[16*tree->mod->ns+13] =    3 ;
834      tree->step_mat[16*tree->mod->ns+14] =    2 ;
835      tree->step_mat[16*tree->mod->ns+15] =    2 ;
836      tree->step_mat[16*tree->mod->ns+16] =    0 ;
837      tree->step_mat[16*tree->mod->ns+17] =    3 ;
838      tree->step_mat[16*tree->mod->ns+18] =    3 ;
839      tree->step_mat[16*tree->mod->ns+19] =    3 ;
840      tree->step_mat[17*tree->mod->ns+ 0] =    3 ;
841      tree->step_mat[17*tree->mod->ns+ 1] =    2 ;
842      tree->step_mat[17*tree->mod->ns+ 2] =    3 ;
843      tree->step_mat[17*tree->mod->ns+ 3] =    3 ;
844      tree->step_mat[17*tree->mod->ns+ 4] =    1 ;
845      tree->step_mat[17*tree->mod->ns+ 5] =    2 ;
846      tree->step_mat[17*tree->mod->ns+ 6] =    2 ;
847      tree->step_mat[17*tree->mod->ns+ 7] =    2 ;
848      tree->step_mat[17*tree->mod->ns+ 8] =    3 ;
849      tree->step_mat[17*tree->mod->ns+ 9] =    3 ;
850      tree->step_mat[17*tree->mod->ns+10] =    2 ;
851      tree->step_mat[17*tree->mod->ns+11] =    2 ;
852      tree->step_mat[17*tree->mod->ns+12] =    2 ;
853      tree->step_mat[17*tree->mod->ns+13] =    2 ;
854      tree->step_mat[17*tree->mod->ns+14] =    3 ;
855      tree->step_mat[17*tree->mod->ns+15] =    2 ;
856      tree->step_mat[17*tree->mod->ns+16] =    3 ;
857      tree->step_mat[17*tree->mod->ns+17] =    0 ;
858      tree->step_mat[17*tree->mod->ns+18] =    2 ;
859      tree->step_mat[17*tree->mod->ns+19] =    3 ;
860      tree->step_mat[18*tree->mod->ns+ 0] =    3 ;
861      tree->step_mat[18*tree->mod->ns+ 1] =    3 ;
862      tree->step_mat[18*tree->mod->ns+ 2] =    2 ;
863      tree->step_mat[18*tree->mod->ns+ 3] =    2 ;
864      tree->step_mat[18*tree->mod->ns+ 4] =    2 ;
865      tree->step_mat[18*tree->mod->ns+ 5] =    2 ;
866      tree->step_mat[18*tree->mod->ns+ 6] =    2 ;
867      tree->step_mat[18*tree->mod->ns+ 7] =    3 ;
868      tree->step_mat[18*tree->mod->ns+ 8] =    2 ;
869      tree->step_mat[18*tree->mod->ns+ 9] =    3 ;
870      tree->step_mat[18*tree->mod->ns+10] =    2 ;
871      tree->step_mat[18*tree->mod->ns+11] =    2 ;
872      tree->step_mat[18*tree->mod->ns+12] =    3 ;
873      tree->step_mat[18*tree->mod->ns+13] =    2 ;
874      tree->step_mat[18*tree->mod->ns+14] =    3 ;
875      tree->step_mat[18*tree->mod->ns+15] =    2 ;
876      tree->step_mat[18*tree->mod->ns+16] =    3 ;
877      tree->step_mat[18*tree->mod->ns+17] =    2 ;
878      tree->step_mat[18*tree->mod->ns+18] =    0 ;
879      tree->step_mat[18*tree->mod->ns+19] =    3 ;
880      tree->step_mat[19*tree->mod->ns+ 0] =    2 ;
881      tree->step_mat[19*tree->mod->ns+ 1] =    3 ;
882      tree->step_mat[19*tree->mod->ns+ 2] =    3 ;
883      tree->step_mat[19*tree->mod->ns+ 3] =    2 ;
884      tree->step_mat[19*tree->mod->ns+ 4] =    3 ;
885      tree->step_mat[19*tree->mod->ns+ 5] =    3 ;
886      tree->step_mat[19*tree->mod->ns+ 6] =    2 ;
887      tree->step_mat[19*tree->mod->ns+ 7] =    2 ;
888      tree->step_mat[19*tree->mod->ns+ 8] =    3 ;
889      tree->step_mat[19*tree->mod->ns+ 9] =    2 ;
890      tree->step_mat[19*tree->mod->ns+10] =    2 ;
891      tree->step_mat[19*tree->mod->ns+11] =    3 ;
892      tree->step_mat[19*tree->mod->ns+12] =    2 ;
893      tree->step_mat[19*tree->mod->ns+13] =    2 ;
894      tree->step_mat[19*tree->mod->ns+14] =    3 ;
895      tree->step_mat[19*tree->mod->ns+15] =    3 ;
896      tree->step_mat[19*tree->mod->ns+16] =    3 ;
897      tree->step_mat[19*tree->mod->ns+17] =    3 ;
898      tree->step_mat[19*tree->mod->ns+18] =    3 ;
899      tree->step_mat[19*tree->mod->ns+19] =    0 ;
900    }
901  else
902    {
903      tree->step_mat[0*tree->mod->ns+0] = 0;
904      tree->step_mat[0*tree->mod->ns+1] = 2;
905      tree->step_mat[0*tree->mod->ns+2] = 1;
906      tree->step_mat[0*tree->mod->ns+3] = 2;
907
908      tree->step_mat[1*tree->mod->ns+0] = 2;
909      tree->step_mat[1*tree->mod->ns+1] = 0;
910      tree->step_mat[1*tree->mod->ns+2] = 2;
911      tree->step_mat[1*tree->mod->ns+3] = 1;
912
913      tree->step_mat[2*tree->mod->ns+0] = 1;
914      tree->step_mat[2*tree->mod->ns+1] = 2;
915      tree->step_mat[2*tree->mod->ns+2] = 0;
916      tree->step_mat[2*tree->mod->ns+3] = 2;
917
918      tree->step_mat[3*tree->mod->ns+0] = 2;
919      tree->step_mat[3*tree->mod->ns+1] = 1;
920      tree->step_mat[3*tree->mod->ns+2] = 2;
921      tree->step_mat[3*tree->mod->ns+3] = 0;
922
923/*       tree->step_mat[0*tree->mod->ns+0] = 0; */
924/*       tree->step_mat[0*tree->mod->ns+1] = 1; */
925/*       tree->step_mat[0*tree->mod->ns+2] = 1; */
926/*       tree->step_mat[0*tree->mod->ns+3] = 1; */
927
928/*       tree->step_mat[1*tree->mod->ns+0] = 1; */
929/*       tree->step_mat[1*tree->mod->ns+1] = 0; */
930/*       tree->step_mat[1*tree->mod->ns+2] = 1; */
931/*       tree->step_mat[1*tree->mod->ns+3] = 1; */
932
933/*       tree->step_mat[2*tree->mod->ns+0] = 1; */
934/*       tree->step_mat[2*tree->mod->ns+1] = 1; */
935/*       tree->step_mat[2*tree->mod->ns+2] = 0; */
936/*       tree->step_mat[2*tree->mod->ns+3] = 1; */
937
938/*       tree->step_mat[3*tree->mod->ns+0] = 1; */
939/*       tree->step_mat[3*tree->mod->ns+1] = 1; */
940/*       tree->step_mat[3*tree->mod->ns+2] = 1; */
941/*       tree->step_mat[3*tree->mod->ns+3] = 0; */
942
943    }
944 
945  For(i,tree->mod->ns) tree->step_mat[i*tree->mod->ns+i] = 0;
946}
947
948/*********************************************************/
949
950/*********************************************************/
Note: See TracBrowser for help on using the repository browser.