website/archive/binaries/mac/src/tcoffee/t_coffee_source/util_declare.c

   1 #include <stdio.h>
   2 #include <stdlib.h>
   3 #include <strings.h>
   4 #include <string.h>
   5 #include <limits.h>
   6 #include <float.h>
   7 #include <assert.h>
   8 #include <math.h>
   9 #include <stdarg.h>
  10
  11 #include "io_lib_header.h"
  12 #include "util_lib_header.h"
  13 #include "define_header.h"
  14 #include "dp_lib_header.h"
  15
  16 void free_pair_wise()
  17 {
  18   //Free static allocated memory
  19   free_proba_pair_wise();
  20 }
  21
  22 /************************************************************************/
  23 /*                                                                      */
  24 /*            CONSTRAINT_LIST                                           */
  25 /*                                                                      */
  26 /*                                                                      */
  27 /************************************************************************/
  28 int *** duplicate_residue_index (int ***r)
  29 {
  30   int a,b,c;
  31   int d1,d2,d3;
  32
  33   int ***nr;
  34
  35   d1=read_array_size_new(r);
  36   nr=vcalloc ( d1, sizeof (int**));
  37   for (a=0; a<d1; a++)
  38     {
  39       d2=read_array_size_new (r[a])-1;
  40       nr[a]=vcalloc ( d2+1, sizeof (int*));
  41       for (b=0; b<d2; b++)
  42         {
  43           d3=read_array_size_new (r[a][b]);
  44           for (c=0; c<d3; c++)nr[a][b][c]=r[a][b][c];
  45         }
  46     }
  47   return nr;
  48 }
  49 int *** declare_residue_index (Sequence *S)
  50 {
  51   int ***r;
  52   int a,b,c;
  53
  54   if ( !S)return NULL;
  55   r=vcalloc ( S->nseq, sizeof (int**));
  56   for ( a=0; a<S->nseq; a++)
  57     {
  58       r[a]=vcalloc (S->len[a]+2, sizeof (int*));//The empty terminator makes a scan possible without knowing len
  59       for ( b=0; b<=S->len[a]; b++)
  60         {
  61           r[a][b]=vcalloc ( 1, sizeof (int));
  62           r[a][b][0]=1;
  63         }
  64     }
  65   return r;
  66 }
  67
  68 Constraint_list * declare_constraint_list_simple ( Sequence *S)
  69 {
  70   return declare_constraint_list (S, NULL, NULL, 0, NULL, NULL);
  71 }
  72
  73 Constraint_list * declare_constraint_list ( Sequence *S, char *name, int *L, int ne,FILE *fp, int **M)
  74     {
  75     Constraint_list *CL;
  76
  77     CL=vcalloc (1, sizeof ( Constraint_list));
  78
  79
  80     CL->S=S;
  81     CL->M=M;
  82
  83     if ( name!=NULL)
  84         {
  85         sprintf ( CL->list_name, "%s", name);
  86
  87         }
  88     CL->cpu=1;
  89     CL->fp=fp;
  90     if (L)
  91       {
  92         HERE ("The USE of L is now Deprecated with Constraint Lists");
  93         exit (0);
  94       }
  95     CL->ne=ne;
  96     CL->entry_len=LIST_N_FIELDS;
  97     CL->el_size=sizeof (CLIST_TYPE);
  98     CL->matrices_list=declare_char(20,20);
  99
 100
 101     CL->weight_field=WE;
 102     if ( S)CL->seq_for_quadruplet=vcalloc ( S->nseq, sizeof (int));
 103     CL->Prot_Blast=vcalloc ( 1, sizeof ( Blast_param));
 104     CL->DNA_Blast=vcalloc ( 1, sizeof ( Blast_param));
 105     CL->Pdb_Blast=vcalloc ( 1, sizeof ( Blast_param));
 106     CL->TC=vcalloc (1, sizeof (TC_param));
 107
 108     //New data structure
 109     CL->residue_index=declare_residue_index (S);
 110
 111
 112
 113     return CL;
 114     }
 115
 116 Constraint_list *free_constraint_list4lib_computation (Constraint_list *CL)
 117 {
 118   if (!CL)return NULL;
 119
 120   free_arrayN(CL->residue_index, 3);
 121   free_int (CL->M, -1);
 122
 123   vfree (CL);
 124   return CL;
 125 }
 126 Constraint_list *duplicate_constraint_list4lib_computation (Constraint_list *CL)
 127 {
 128   Constraint_list *SCL;
 129   SCL=vcalloc (1, sizeof ( Constraint_list));
 130   SCL[0]=CL[0];
 131   SCL->S=CL->S;
 132   SCL->RunName=CL->RunName;
 133
 134   SCL->max_L_len=0;
 135   SCL->M=NULL;
 136   SCL->ne=0;
 137
 138   return SCL;
 139 }
 140 Constraint_list *duplicate_constraint_list_soft (Constraint_list *CL)
 141 {
 142   /*Duplication that does not copy the long lists*/
 143   return copy_constraint_list (CL,SOFT_COPY);
 144 }
 145 Constraint_list *duplicate_constraint_list (Constraint_list *CL)
 146      {
 147      /*Duplicate everything in the constraint_list*/
 148      return copy_constraint_list (CL,HARD_COPY);
 149      }
 150 Constraint_list *copy_constraint_list (Constraint_list *CL, int mode)
 151     {
 152     Constraint_list *NCL;
 153     Sequence *S;
 154     int a, b;
 155
 156
 157
 158     /*Sequences*/
 159
 160
 161       S=(mode==HARD_COPY)?duplicate_sequence (CL->S):CL->S;
 162
 163
 164       if (mode==HARD_COPY)
 165         NCL=declare_constraint_list (S, NULL, NULL,0, NULL, NULL);
 166       else
 167         {
 168           NCL=vcalloc ( 1, sizeof (Constraint_list));
 169           NCL[0]=CL[0];
 170         }
 171
 172
 173       NCL->copy_mode=mode;
 174       if (mode==SOFT_COPY)NCL->pCL=CL;
 175       NCL->S=S;
 176       /*master*/
 177       if (mode==HARD_COPY && CL->master)
 178         {NCL->master=vcalloc ( S->nseq, sizeof(int));
 179         for ( a=0; a< S->nseq; a++)
 180           NCL->master[a]=CL->master[a];
 181         }
 182       else if (mode==SOFT_COPY)
 183         {
 184           NCL->seq_for_quadruplet=CL->seq_for_quadruplet;
 185         }
 186       NCL->o2a_byte=CL->o2a_byte;
 187
 188       /*struc List*/
 189       NCL->STRUC_LIST=(mode==HARD_COPY)?duplicate_sequence (CL->STRUC_LIST):CL->STRUC_LIST;
 190       sprintf ( NCL->align_pdb_param_file, "%s", CL->align_pdb_param_file);
 191       sprintf ( NCL->align_pdb_hasch_mode, "%s", CL->align_pdb_hasch_mode);
 192
 193
 194       NCL->W=(mode==HARD_COPY)?duplicate_weights (CL->W):CL->W;
 195       NCL->DM=(mode==HARD_COPY)?duplicate_distance_matrix (CL->DM):CL->DM;
 196       NCL->ktupDM=(mode==HARD_COPY)?duplicate_distance_matrix (CL->ktupDM):CL->ktupDM;
 197       NCL->RunName=CL->RunName;
 198
 199       if (  mode==HARD_COPY && CL->translation){NCL->translation=vcalloc ((CL->S)->nseq, sizeof (int)); for ( a=0; a< (CL->S)->nseq; a++)NCL->translation[a]=CL->translation[a];}
 200       else{NCL->translation=CL->translation;}
 201
 202       NCL->out_aln_format=(mode==HARD_COPY)?duplicate_char (CL->out_aln_format, -1, -1):CL->out_aln_format;
 203       NCL->n_out_aln_format=CL->n_out_aln_format;
 204
 205     /*Packing Sequence: To use with domain analysis*/
 206       NCL->packed_seq_lu=(mode==HARD_COPY)?duplicate_int (CL->packed_seq_lu, -1, -1):CL->packed_seq_lu;
 207     /*DATA*/
 208       if (CL->fp)(mode==HARD_COPY)?NCL->fp=vtmpfile():CL->fp;
 209
 210       if ( mode==HARD_COPY)
 211         {
 212           NCL->residue_index=duplicate_residue_index (NCL->residue_index);
 213         }
 214       else NCL->residue_index=CL->residue_index;
 215
 216
 217      if ( mode==HARD_COPY)
 218        {
 219          NCL->M=copy_int ( CL->M,NCL->M,-1, -1);
 220        }
 221      else
 222        NCL->M=CL->M;
 223
 224
 225     /*List Information*/
 226       NCL->ne=CL->ne;
 227       sprintf ( NCL->list_name, "%s", CL->list_name);
 228       NCL->entry_len=CL->entry_len;
 229       NCL->el_size=CL->el_size;
 230
 231     /*Normalisation information*/
 232       NCL->filter_lib=CL->filter_lib;
 233       NCL->normalise=CL->normalise;
 234       NCL->overweight=CL->overweight;
 235       NCL->max_ext_value=CL->max_ext_value;
 236       NCL->max_value=CL->max_value;
 237
 238     /*Pair wise alignment method*/
 239       NCL->pw_parameters_set=CL->pw_parameters_set;
 240       NCL->gop=CL->gop;
 241       NCL->f_gop=CL->f_gop;
 242       NCL->gep=CL->gep;
 243       NCL->f_gep=CL->f_gep;
 244
 245       NCL->nomatch=CL->nomatch;
 246
 247       NCL->TG_MODE=CL->TG_MODE;
 248       NCL->F_TG_MODE=CL->F_TG_MODE;
 249
 250       sprintf ( NCL->dp_mode, "%s", CL->dp_mode);
 251       NCL->maximise=CL->maximise;
 252       sprintf ( NCL->matrix_for_aa_group, "%s", CL->matrix_for_aa_group);
 253       sprintf ( NCL->method_matrix, "%s", CL->method_matrix);
 254
 255       NCL->diagonal_threshold=CL->diagonal_threshold;
 256       NCL->ktup=CL->ktup;
 257
 258       NCL->use_fragments=CL->use_fragments;
 259       NCL->fasta_step=CL->fasta_step;
 260       NCL->lalign_n_top=CL->lalign_n_top;
 261       NCL->sw_min_dist=CL->sw_min_dist;
 262       NCL->matrices_list=(mode==HARD_COPY)?duplicate_char (CL->matrices_list, -1, -1):CL->matrices_list;
 263       NCL->n_matrices=CL->n_matrices;
 264
 265       sprintf (NCL->distance_matrix_mode, "%s", CL->distance_matrix_mode);
 266       sprintf (NCL->distance_matrix_sim_mode, "%s", CL->distance_matrix_sim_mode);
 267
 268       sprintf (NCL->tree_mode, "%s", CL->tree_mode);
 269       NCL->tree_aln=(mode==HARD_COPY)?copy_aln (CL->tree_aln, NULL):CL->tree_aln;
 270     /*Functions used for dynamic programming and Evaluation*/
 271       NCL->no_overaln=CL->no_overaln;
 272       NCL->profile_mode=CL->profile_mode;
 273       sprintf ( NCL->profile_comparison, "%s",CL->profile_comparison);
 274       NCL->get_dp_cost=CL->get_dp_cost;
 275       NCL->evaluate_residue_pair=CL->evaluate_residue_pair;
 276       NCL->pair_wise=CL->pair_wise;
 277
 278       NCL->weight_field=CL->weight_field;
 279       NCL->max_n_pair=CL->max_n_pair;
 280
 281     /*threading parameters*/
 282       NCL->Prot_Blast=(mode==HARD_COPY)?duplicate_blast_param ( CL->Prot_Blast):CL->Prot_Blast;
 283       NCL->DNA_Blast =(mode==HARD_COPY)?duplicate_blast_param ( CL->DNA_Blast):CL->DNA_Blast;
 284       NCL->Pdb_Blast =(mode==HARD_COPY)?duplicate_blast_param ( CL->Pdb_Blast):CL->Pdb_Blast;
 285       NCL->TC =(mode==HARD_COPY)?duplicate_TC_param ( CL->TC):CL->TC;
 286
 287     /*Split parameters*/
 288       NCL->split=CL->split;
 289       NCL->split_nseq_thres= CL->split_nseq_thres;
 290       NCL->split_score_thres= CL->split_score_thres;
 291     /*Structural status*/
 292       NCL->check_pdb_status=CL->check_pdb_status;
 293     /*log*/
 294       sprintf ( NCL->method_log, "%s",CL->method_log);
 295       sprintf ( NCL->evaluate_mode, "%s",CL->evaluate_mode);
 296     /* Gene Prediction*/
 297       sprintf ( NCL->genepred_score, "%s",CL->genepred_score);
 298
 299     /*Parameters for domain extraction*/
 300       NCL->moca=(mode==HARD_COPY)?duplicate_moca ( CL->moca):CL->moca;
 301
 302
 303
 304     /*Functions for hiding forbiden pairs of residues*/
 305       /* Copy only for soft_copy*/
 306       if (mode==SOFT_COPY)
 307         {
 308           NCL->forbiden_pair_list=CL->forbiden_pair_list;
 309         }
 310    /*extention properties:*/
 311       NCL->nseq_for_quadruplet=CL->nseq_for_quadruplet;
 312       if (mode==HARD_COPY && CL->seq_for_quadruplet)
 313         {NCL->seq_for_quadruplet=vcalloc ( S->nseq, sizeof(int));
 314         for ( a=0; a< S->nseq; a++)
 315           NCL->seq_for_quadruplet[a]=CL->seq_for_quadruplet[a];
 316         }
 317       else if (mode==SOFT_COPY)
 318         {
 319           NCL->seq_for_quadruplet=CL->seq_for_quadruplet;
 320         }
 321
 322    /*extention properties: Do only a soft copy*/
 323       /* Not To be copied yet */
 324       if ( mode==SOFT_COPY)
 325         {
 326           NCL->extend_jit=CL->extend_jit;
 327           NCL->extend_threshold=CL->extend_threshold;
 328           sprintf ( NCL->extend_clean_mode, "%s", CL->extend_clean_mode);
 329           sprintf ( NCL->extend_compact_mode, "%s", CL->extend_compact_mode);
 330         }
 331
 332     /*Lookup table parameteres*/
 333       NCL->chunk= CL->chunk;
 334       /* Do NOT copy NCL->seq_indexed, NCL->start_index, NCL->max_L_len, NCL->chunk*/
 335       /*
 336         if ( mode==SOFT_COPY)
 337         {
 338           NCL->seq_indexed=CL->seq_indexed;
 339           NCL->start_index=CL->start_index;
 340           NCL->end_index=CL->start_index;
 341           NCL->max_L_len=CL->max_L_len;
 342           }
 343       */
 344     /*PDB STRUCTURE ALIGNMENTS*/
 345       /* Do only a soft copy */
 346       if ( mode==SOFT_COPY)
 347         {
 348           NCL->T=CL->T;
 349         }
 350     /*MISC*/
 351        NCL->cpu=CL->cpu;
 352        NCL->local_stderr=CL->local_stderr;
 353        sprintf (NCL->multi_thread, "%s", CL->multi_thread);
 354
 355     return NCL;
 356     }
 357 Constraint_list *free_constraint_list_full (Constraint_list *CL)
 358 {
 359   free_sequence (free_constraint_list (CL), -1);
 360   return NULL;
 361 }
 362 Sequence *free_constraint_list (Constraint_list *CL)
 363     {
 364     Sequence *S;
 365     int a, b;
 366     Constraint_list *pCL;
 367
 368
 369     /*Prepare the selective freeing of the CL data structure:
 370       If the CL has been obtained from copy, every pointer that is identical to the parent CL (CL->pCL)
 371       will not be saved.
 372     */
 373
 374
 375     if ( !CL)return NULL;
 376     else S=CL->S;
 377
 378     if ( CL->copy_mode==SOFT_COPY && !CL->pCL)
 379       {vfree(CL); return S;}
 380     else if ( CL->copy_mode==SOFT_COPY)
 381       {
 382
 383         pCL=CL->pCL;
 384         CL->residue_index=NULL;
 385
 386         if ( CL->M                      ==pCL->M                       )CL->M=NULL;
 387
 388         if (CL->start_index             ==pCL->start_index             )CL->start_index=NULL;
 389         if (CL->end_index             ==pCL->end_index                 )CL->end_index=NULL;
 390
 391         if ( CL->fp                     ==pCL->fp                      )CL->fp=NULL;
 392         if ( CL->matrices_list          ==pCL->matrices_list           )CL->matrices_list=NULL;
 393
 394
 395         if ( CL->STRUC_LIST             ==pCL->STRUC_LIST              )CL->STRUC_LIST=NULL;
 396         if ( CL->W                      ==pCL->W                       )CL->W=NULL;
 397         if ( CL->DM                     ==pCL->DM                      )CL->DM=NULL;
 398         if ( CL->ktupDM                 ==pCL->ktupDM                      )CL->ktupDM=NULL;
 399
 400
 401         if ( CL->translation            ==pCL->translation             )CL->translation=NULL;
 402         if ( CL->moca                   ==pCL->moca                    )CL->moca=NULL;
 403         if ( CL->Prot_Blast             ==pCL->Prot_Blast              )CL->Prot_Blast=NULL;
 404         if ( CL->DNA_Blast              ==pCL->DNA_Blast               )CL->DNA_Blast=NULL;
 405         if ( CL->Pdb_Blast              ==pCL->Pdb_Blast               )CL->Pdb_Blast=NULL;
 406         if ( CL->seq_for_quadruplet     ==pCL->seq_for_quadruplet      )CL->seq_for_quadruplet=NULL;
 407         if ( CL->TC                      ==pCL->TC                       )CL->TC=NULL;
 408
 409       }
 410
 411
 412     /*End of selective freeing of the CL data structure*/
 413
 414
 415
 416     if ( CL->residue_index)free_arrayN(CL->residue_index, 3);
 417
 418     if ( CL->M)free_int (CL->M, -1);
 419     if ( CL->fp)vfclose (CL->fp);
 420     if ( CL->matrices_list)free_char(CL->matrices_list,-1);
 421
 422
 423     if ( CL->start_index)free_int ( CL->start_index,-1);
 424     if ( CL->end_index)free_int ( CL->end_index,-1);
 425
 426
 427
 428
 429     if ( CL->STRUC_LIST)free_sequence ( CL->STRUC_LIST, (CL->STRUC_LIST)->nseq);
 430     if ( CL->W)free_weights (CL->W);
 431
 432     CL->DM=free_distance_matrix (CL->DM);
 433     CL->ktupDM=free_distance_matrix (CL->ktupDM);
 434
 435     if ( CL->translation)vfree(CL->translation);
 436     if ( CL->moca)free_moca (CL->moca);
 437     if ( CL->Prot_Blast)free_blast_param ( CL->Prot_Blast);
 438     if ( CL->DNA_Blast) free_blast_param ( CL->DNA_Blast);
 439     if ( CL->Pdb_Blast) free_blast_param ( CL->Pdb_Blast);
 440     if ( CL->TC) free_TC_param ( CL->TC);
 441
 442     if (CL->seq_for_quadruplet)vfree (CL->seq_for_quadruplet);
 443
 444     vfree(CL);
 445     return S;
 446     }
 447
 448 Distance_matrix * free_distance_matrix ( Distance_matrix *DM)
 449 {
 450   if (!DM)return NULL;
 451   free_int ( DM->similarity_matrix,-1);
 452   free_int ( DM->distance_matrix,-1);
 453   free_int ( DM->score_similarity_matrix,-1);
 454   vfree (DM);
 455   return NULL;
 456 }
 457 Distance_matrix * duplicate_distance_matrix ( Distance_matrix *DMin)
 458 {
 459   Distance_matrix *DM;
 460   if (!DMin) return NULL;
 461
 462   DM=vcalloc ( 1, sizeof (Distance_matrix));
 463   DM->similarity_matrix=duplicate_int ( DMin->similarity_matrix, -1, -1);
 464   DM->distance_matrix=duplicate_int ( DMin->distance_matrix, -1, -1);
 465   DM->score_similarity_matrix=duplicate_int ( DMin->score_similarity_matrix, -1, -1);
 466   return DM;
 467 }
 468
 469 /************************************************************************/
 470 /*                                                                      */
 471 /*            MOCA Functions                                            */
 472 /*                                                                      */
 473 /*                                                                      */
 474 /************************************************************************/
 475 Moca * duplicate_moca ( Moca *m)
 476       {
 477         Moca *nm;
 478
 479         if ( m==NULL)return m;
 480
 481         nm=vcalloc ( 1, sizeof (Moca));
 482
 483         nm->moca_scale=m->moca_scale;
 484         nm->evaluate_domain=m->evaluate_domain;
 485         nm->moca_threshold=m->moca_threshold;
 486         nm->cache_cl_with_domain=m->cache_cl_with_domain;
 487         if ( m->forbiden_residues)nm->forbiden_residues=copy_int  (m->forbiden_residues,nm->forbiden_residues,    -1, -1);
 488         nm->make_nol_aln=m->make_nol_aln;
 489
 490
 491         return nm;
 492       }
 493 Moca * free_moca ( Moca *m)
 494       {
 495         if ( m->forbiden_residues)free_int ( m->forbiden_residues, -1);
 496         vfree ( m);
 497         return NULL;
 498       }
 499 /************************************************************************/
 500 /*                                                                      */
 501 /*            TC_param Functions                                            */
 502 /*                                                                      */
 503 /*                                                                      */
 504 /************************************************************************/
 505 TC_param * duplicate_TC_param ( TC_param*B)
 506 {
 507   TC_param *N;
 508   N=vcalloc (1, sizeof ( TC_param));
 509   memcpy(B, N, sizeof(TC_param));
 510   return N;
 511   }
 512 TC_param * free_TC_param ( TC_param*B)
 513 {
 514   vfree (B);
 515   return NULL;
 516 }
 517 /************************************************************************/
 518 /*                                                                      */
 519 /*            Blast_param Functions                                            */
 520 /*                                                                      */
 521 /*                                                                      */
 522 /************************************************************************/
 523 Blast_param * duplicate_blast_param ( Blast_param*B)
 524 {
 525   Blast_param *N;
 526   N=vcalloc (1, sizeof ( Blast_param));
 527   sprintf ( N->blast_server, "%s", B->blast_server);
 528   sprintf ( N->db, "%s", B->db);
 529   N->min_id=B->min_id;
 530   N->max_id=B->min_id;
 531   N->min_cov=B->min_cov;
 532   return N;
 533 }
 534 Blast_param * free_blast_param ( Blast_param*B)
 535 {
 536   vfree (B);
 537   return NULL;
 538 }
 539
 540 /************************************************************************/
 541 /*                                                                      */
 542 /*            PDB Functions                                             */
 543 /*                                                                      */
 544 /*                                                                      */
 545 /************************************************************************/
 546 Structure* declare_structure ( int n, char **array)
 547     {
 548     Structure *S;
 549     int a;
 550
 551     S=vcalloc (1, sizeof (Structure));
 552     S->n_fields=1;
 553     S->nseq=n;
 554
 555     S->struc=vcalloc ( n, sizeof (int**));
 556     S->len=vcalloc ( n, sizeof (int));
 557     for ( a=0; a< n; a++)
 558         {
 559         S->len[a]=strlen(array[a]);
 560         S->struc[a]=declare_int ( strlen ( array[a])+2, 1);
 561         }
 562     return S;
 563     }
 564
 565 Structure *extend_structure ( Structure *S)
 566     {
 567     int a, b;
 568
 569
 570     for ( a=0; a< S->nseq; a++)
 571         {
 572         for ( b=0; b< S->len[a]; b++)
 573                 S->struc[a][b]=vrealloc ( S->struc[a][b],( S->n_fields+1)*sizeof (int));
 574         }
 575     S->n_fields++;
 576     return S;
 577     }
 578 Sequence * reset_sequence_len (Sequence *S)
 579 {
 580   int min,max,a,l;
 581
 582   if ( !S || !S->nseq)return 0;
 583   max=min=strlen (S->seq[0]);
 584   for (a=0; a<S->nseq; a++)
 585     {
 586       l=strlen (S->seq[a]);
 587       min=MIN(l,min);
 588       max=MAX(l,max);
 589       S->len[a]=l;
 590     }
 591   S->max_len=max;
 592   S->min_len=min;
 593   return S;
 594 }
 595 Sequence * declare_sequence ( int min, int max, int nseq)
 596     {
 597     Sequence *LS;
 598
 599
 600
 601     LS=vcalloc (1, sizeof ( Sequence));
 602
 603     LS->seq_comment=declare_char ( nseq,COMMENT_SIZE);
 604     LS->aln_comment=declare_char ( nseq,COMMENT_SIZE);
 605
 606     LS->file=declare_char( nseq,STRING+1);
 607     LS->seq=declare_char ( nseq, max+1);
 608     LS->name=declare_char( nseq,MAXNAMES+1);
 609
 610     LS->len=vcalloc ( nseq, sizeof (int));
 611     LS->max_len=max;
 612     LS->min_len=min;
 613     LS->nseq=nseq;
 614     LS->max_nseq=nseq;
 615     LS->type=vcalloc(30, sizeof (char));
 616     LS->T=declare_arrayN(2, sizeof (Template), nseq, 1);
 617
 618
 619     LS->dc=declare_int (nseq, 2);
 620     return LS;
 621     }
 622 Sequence * realloc_sequence   (Sequence *OUT, int new_nseq, int max_len)
 623     {
 624
 625
 626       if ( new_nseq<OUT->max_nseq)return OUT;
 627
 628       OUT->min_len =MIN(OUT->min_len,max_len);
 629       OUT->max_len =MAX(OUT->max_len,max_len);
 630       OUT->seq_comment =new_realloc_char ( OUT->seq_comment, new_nseq,COMMENT_SIZE);
 631       OUT->aln_comment =new_realloc_char ( OUT->aln_comment, new_nseq,COMMENT_SIZE);
 632
 633       OUT->seq     =new_realloc_char ( OUT->seq,     new_nseq,OUT->max_len+1);
 634       OUT->name    =new_realloc_char ( OUT->name,    new_nseq,MAXNAMES+1);
 635
 636
 637       OUT->file    =new_realloc_char ( OUT->file,    new_nseq,STRING+1);
 638       OUT->len     =vrealloc     ( OUT->len,     (new_nseq+1)*sizeof (int));
 639
 640       OUT->T=(Template**)realloc_arrayN (2, (void **)OUT->T,sizeof (Template), new_nseq, 1);
 641       OUT->dc=(int **)realloc_arrayN (2, (void **)OUT->dc,sizeof (int), new_nseq, 2);
 642
 643       OUT->max_nseq=new_nseq;
 644       return OUT;
 645     }
 646
 647 Sequence * duplicate_sequence (Sequence *S )
 648     {
 649     Sequence *LS;
 650     int a, b;
 651
 652
 653     if (S==NULL)return S;
 654     LS=declare_sequence (S->min_len, S->max_len, S->nseq);
 655     for (b=0, a=0; a<S->nseq; a++)
 656         {
 657           if (S->seq && S->seq[a])
 658             {
 659
 660               sprintf ( LS->file[b], "%s", S->file[a]);
 661               if ( S->seq_comment && S->seq_comment[a])sprintf ( LS->seq_comment[b], "%s", S->seq_comment[a]);
 662               if ( S->aln_comment && S->aln_comment[a])sprintf ( LS->aln_comment[b], "%s", S->aln_comment[a]);
 663               if ( S->seq && S->seq[a])sprintf  ( LS->seq[b],  "%s", S->seq[a]);
 664               if ( S->name&& S->name[a])sprintf ( LS->name[b], "%s", S->name[a]);
 665               LS->dc[b][0]=S->dc[a][0];
 666               LS->dc[b][1]=S->dc[a][1];
 667               LS->len[b]=S->len[a];
 668               LS->T[b][0]=S->T[a][0];
 669
 670               b++;
 671
 672             }
 673         }
 674
 675     LS->max_len=S->max_len;
 676     LS->min_len=S->min_len;
 677     LS->nseq=b;
 678
 679     if (S->W)LS->W=duplicate_weights (S->W);
 680     sprintf ( LS->type, "%s", S->type);
 681     sprintf ( LS->template_file, "%s", S->template_file);
 682     LS->max_nseq=S->nseq;
 683
 684     return LS;
 685     }
 686
 687 void free_sequence ( Sequence *LS, int nseq)
 688         {
 689
 690
 691         if ( !LS) return;
 692
 693
 694         free_char ( LS->file, -1);
 695         free_char ( LS->seq_comment, -1);
 696         free_char ( LS->aln_comment, -1);
 697         free_char ( LS->seq, -1);
 698         free_char ( LS->name,-1);
 699         free_int  (LS->dc, -1);
 700         free_arrayN((void*)LS->T, 2);
 701         vfree (LS->type);
 702         vfree (LS->len);
 703         free_weights (LS->W);
 704
 705         vfree (LS);
 706
 707         }
 708 /************************************************************************/
 709 /*                                                                      */
 710 /*            Weights Functions                                         */
 711 /*                                                                      */
 712 /*                                                                      */
 713 /************************************************************************/
 714 Weights* declare_weights ( int nseq)
 715         {
 716         Weights *W;
 717
 718         W=vcalloc ( 1, sizeof ( Weights));
 719         W->comments=vcalloc ( 1000, sizeof (char));
 720         W->nseq=nseq;
 721         W->mode=vcalloc (FILENAMELEN, sizeof (char));
 722         W->seq_name= declare_char ( W->nseq*2, 200);
 723         W->PW_SD=declare_float ( W->nseq, W->nseq);
 724         W->PW_ID=declare_float ( W->nseq, W->nseq);
 725         W->SEQ_W=vcalloc ( W->nseq, sizeof ( float));
 726         return W;
 727         }
 728 Weights* duplicate_weights (Weights *W)
 729      {
 730        Weights *NW;
 731        int a, b, c;
 732
 733        NW=declare_weights (W->nseq);
 734        sprintf ( NW->comments, "%s", W->comments);
 735        sprintf ( NW->mode, "%s", W->mode);
 736        for (a=0, c=0; a< W->nseq; a++)
 737          {
 738            if ( W->seq_name[a])
 739              {
 740                sprintf ( NW->seq_name[c], "%s", W->seq_name[a]);
 741                NW->SEQ_W[c]=W->SEQ_W[a];
 742                for(b=0; b< W->nseq; b++)
 743                  {
 744                    NW->PW_SD[c][b]=W->PW_SD[a][b];
 745                    NW->PW_ID[c][b]=W->PW_ID[a][b];
 746                  }
 747                c++;
 748              }
 749          }
 750        return NW;
 751      }
 752 Weights* free_weights ( Weights* W)
 753         {
 754
 755         if ( !W)return NULL;
 756
 757         vfree(W->comments);
 758
 759
 760         vfree(W->mode);
 761
 762         free_char(W->seq_name, -1);
 763         free_float(W->PW_SD,-1);
 764         free_float(W->PW_ID, -1);
 765         vfree(W->SEQ_W);
 766         vfree(W);
 767         return NULL;
 768         }
 769
 770
 771 Alignment* copy_aln ( Alignment *A, Alignment *B)
 772         {
 773           int a, b;
 774           int nnseq;
 775           int nlen;
 776           /*      c[100]=10;*/
 777
 778
 779
 780           if ( A==NULL){free_aln(B); return NULL;}
 781
 782           nnseq=MAX(A->nseq, A->max_n_seq);
 783           nlen=A->len_aln+1;
 784           if (B)
 785             B=realloc_alignment2 (B, nnseq, nlen);
 786           else
 787             B=declare_aln2 (nnseq, nlen);
 788           B->S=A->S;
 789
 790
 791           /*SIZES*/
 792           B->max_len=A->max_len;
 793           B->min_len=A->min_len;
 794           B->declared_len=nlen;
 795           B->max_n_seq=nnseq;
 796
 797           B->nseq=A->nseq;
 798           B->len_aln=A->len_aln;
 799
 800
 801 /*sequence Information*/
 802             if ( A->generic_comment)
 803               {
 804                 vfree(B->generic_comment);
 805                 B->generic_comment=vcalloc (strlen(A->generic_comment)+1, sizeof (char));
 806                 sprintf ( B->generic_comment, "%s", A->generic_comment);
 807               }
 808             if ( (A->S)==NULL){vfree (B->len); B->len=vcalloc ( A->max_n_seq, sizeof (int));}
 809             ga_memcpy_int ( A->len, B->len, B->nseq);
 810
 811             B->seq_comment=copy_char ( A->seq_comment,  B->seq_comment,  -1,-1);
 812             B->aln_comment=copy_char ( A->aln_comment,  B->aln_comment,  -1,-1);
 813
 814             B->name=copy_char ( A->name,     B->name,     -1,-1);
 815
 816             B->file=copy_char ( A->file,     B->file,     -1,-1);
 817             B->tree_order=copy_char ( A->tree_order,     B->tree_order,     -1,-1);
 818             B->expanded_order=A->expanded_order;
 819             free_char ( B->seq_al, -1);
 820             B->seq_al=declare_char(B->max_n_seq, B->declared_len);
 821             // HERE ("A: MAX_NSEQ=%d %d %d %d",B->nseq, B->max_n_seq, B->declared_len, B->len_aln);
 822             // HERE ("B: MAX_NSEQ=%d %d %d %d",A->nseq, A->max_n_seq, A->declared_len, A->len_aln);
 823             for ( a=0; a< nnseq; a++)
 824               {
 825                 if (A->seq_al[a])
 826                   {
 827                     for ( b=0; b< A->len_aln; b++)
 828                       B->seq_al[a][b]=A->seq_al[a][b];
 829                   }
 830               }
 831
 832
 833
 834             B->order=copy_int  ( A->order,    B->order,    -1, -1);
 835             B->S=A->S;
 836             if (A->seq_cache)
 837                 {
 838                 B->seq_cache=copy_int  ( A->seq_cache,    B->seq_cache,-1,-1);
 839                 }
 840
 841             if (A->cdna_cache)
 842                 {
 843                 B->cdna_cache=copy_int  ( A->cdna_cache,    B->cdna_cache,-1,-1);
 844                 }
 845
 846             B->P=copy_profile (A->P);
 847
 848             B->Dp_result=A->Dp_result;
 849
 850 /*Score*/
 851
 852             if ( (A->S)==NULL){vfree (B->score_seq); B->score_seq=vcalloc ( A->max_n_seq, sizeof (int));}
 853             ga_memcpy_int(  A->score_seq,B->score_seq,B->nseq);
 854             B->score_res=A->score_res;
 855
 856             B->score_aln=A->score_aln;
 857             B->score=A->score;
 858             B->ibit=A->ibit;
 859             B->cpu=A->cpu;
 860             B->finished=A->finished;
 861
 862 /*Output Options*/
 863             B->output_res_num=A->output_res_num;
 864             B->residue_case=A->residue_case;
 865             B->expand=A->expand;
 866
 867             B->CL=A->CL;
 868             B->random_tag=A->random_tag;
 869
 870 /*Make the function Recursive */
 871             if ( A->A)
 872               {
 873                 B->A=copy_aln (A->A, NULL);
 874               }
 875             else B->A=NULL;
 876
 877             return B;
 878         }
 879
 880 Alignment* shrink_aln ( Alignment *A, int nseq, int *list)
 881         {
 882         Alignment *B=NULL;
 883         int a,seq;
 884
 885         B=copy_aln (A, B);
 886         for ( a=0; a< nseq; a++)
 887             {
 888             seq=list[a];
 889             sprintf ( A->seq_comment[a], "%s",B->seq_comment[seq]);
 890             sprintf ( A->aln_comment[a], "%s",B->aln_comment[seq]);
 891
 892             sprintf ( A->seq_al [a], "%s",B->seq_al [seq]);
 893             A->order[a][0]=B->order[seq][0];
 894             A->order[a][1]=B->order[seq][1];
 895             A->order[a][2]=B->order[seq][2];
 896             A->order[a][3]=B->order[seq][3];
 897             A->order[a][4]=B->order[seq][4];
 898
 899             A->score_seq[a]=B->score_seq[seq];
 900             A->len[a]=B->len[seq];
 901             }
 902         A->nseq=nseq;
 903         A->len_aln=strlen (A->seq_al[0]);
 904         free_aln (B);
 905         return A;
 906         }
 907 Alignment* extract_sub_aln2 ( Alignment *B, int ns, char **ls)
 908         {
 909           int *list;
 910           Alignment *A;
 911
 912           list=name_array2index_array(ls, ns, B->name, B->nseq);
 913           A=extract_sub_aln ( B,ns, list);
 914           vfree (list);
 915           return A;
 916         }
 917 Alignment* extract_sub_aln ( Alignment *B, int nseq, int *list)
 918         {
 919         Alignment *A=NULL;
 920         int a,b,n,seq;
 921
 922         A=declare_aln2(nseq, B->len_aln+1);
 923         for ( n=0,a=0; a< nseq; a++)
 924             {
 925             seq=list[a];
 926             if ( seq==-1)continue;
 927             else n++;
 928             sprintf ( A->seq_comment[a], "%s",B->seq_comment[seq]);
 929             sprintf ( A->aln_comment[a], "%s",B->aln_comment[seq]);
 930             sprintf ( A->name[a], "%s",B->name[seq]);
 931
 932
 933             for (b=0; b<=B->len_aln; b++)A->seq_al [a][b]=B->seq_al [seq][b];
 934             A->order[a][0]=B->order[seq][0];
 935             A->order[a][1]=B->order[seq][1];
 936             A->order[a][2]=B->order[seq][2];
 937             A->order[a][3]=B->order[seq][3];
 938             A->order[a][4]=B->order[seq][4];
 939
 940             A->score_seq[a]=B->score_seq[seq];
 941             A->len[a]=B->len[seq];
 942             }
 943         A->nseq=n;
 944         A->len_aln=B->len_aln;
 945         return A;
 946         }
 947
 948 Alignment *declare_aln2 ( int nseq, int len)
 949         {
 950           Sequence *S;
 951           Alignment *A;
 952
 953           S=vcalloc ( 1, sizeof ( Sequence));
 954           S->nseq=nseq;
 955           S->max_len=len;
 956
 957           A=declare_aln (S);
 958           A->S=NULL;
 959           vfree(S);
 960           return A;
 961         }
 962
 963
 964
 965 Alignment *declare_aln ( Sequence *S){return declare_Alignment(S);}
 966
 967 Alignment *declare_Alignment ( Sequence *S)
 968         {
 969         Alignment *LA;
 970         int a;
 971
 972         /*ordre:
 973           [x][0]= which is the xth seq of aln
 974           [x][1]= how many deleted residues before the first one
 975         */
 976
 977
 978         LA=vcalloc (1, sizeof ( Alignment));
 979         aln_stack (LA, DECLARE_ALN);
 980         if ( S==NULL)
 981             {
 982               LA->declared_len=MAX_LEN_ALN;
 983               LA->max_n_seq=MAX_N_SEQ;
 984             }
 985         else
 986           {
 987             LA->declared_len=2*S->max_len+1;
 988             LA->max_n_seq=S->nseq+1;
 989           }
 990         LA->S=S;
 991
 992
 993         LA->seq_comment=declare_char (LA->max_n_seq, COMMENT_SIZE);
 994         LA->aln_comment=declare_char (LA->max_n_seq, COMMENT_SIZE);
 995
 996
 997         LA->seq_al=declare_char ( LA->max_n_seq,LA->declared_len );
 998         LA->name=declare_char (LA->max_n_seq, MAXNAMES+1);
 999
1000
1001         LA->file=declare_char (LA->max_n_seq, STRING);
1002         LA->tree_order=declare_char (LA->max_n_seq, STRING);
1003         LA->order= declare_int (LA->max_n_seq , 5);
1004         //order[a][0]: sequence index in S
1005         //order[a][1]: offset of the sequence
1006         //order[a][2]: used by sw_gotoh_pair_wise
1007         //order[a][3]: used by sw_gotoh_pair_wise
1008         //order[a][4]: weight, -1
1009         LA->score_seq= vcalloc (LA->max_n_seq, sizeof (int));
1010
1011         for ( a=0; a< LA->max_n_seq; a++)LA->order[a][0]=a;
1012
1013         LA->len_aln=0;
1014         LA->score_aln=0;
1015         LA->len=vcalloc (LA->max_n_seq, sizeof (int));
1016
1017         if (S && S->name)for ( a=0; a<S->nseq; a++)
1018           {
1019             sprintf ( LA->name[a], "%s", S->name[a]);
1020
1021           }
1022
1023         return LA;
1024
1025         }
1026 Alignment * realloc_aln ( Alignment *A, int new_len){return realloc_alignment(A, new_len);}
1027 Alignment * realloc_alignment ( Alignment *A, int new_len)
1028         {
1029         if (A==NULL)A=declare_Alignment (NULL);
1030
1031         return realloc_alignment2( A, A->max_n_seq,new_len);
1032         }
1033
1034 Alignment * realloc_aln2 ( Alignment *A, int n_nseq, int n_len){return realloc_alignment2(A, n_nseq, n_len);}
1035
1036
1037
1038 Alignment * realloc_alignment2 ( Alignment *A, int n_nseq, int n_len)
1039         {
1040         int a;
1041         int len, nseq;
1042         int delta_len, delta_nseq;
1043
1044         if ( A==NULL) A=declare_Alignment(NULL);
1045
1046         n_len++;
1047         n_nseq++;
1048
1049         len=A->declared_len;
1050         nseq=A->max_n_seq;
1051
1052         n_len=MAX(len, n_len);
1053         n_nseq=MAX(nseq,n_nseq);
1054         delta_nseq=MAX(0,n_nseq-nseq);
1055         delta_len =MAX(0,n_len-len);
1056
1057         if ( delta_nseq<=0 && delta_len<=0)return A;
1058
1059
1060         else
1061             {
1062                 A->len          =vrealloc( A->len      , sizeof (int)*n_nseq);
1063                 for (a=nseq; a< n_nseq; a++)A->len[a]=0;
1064
1065                 A->declared_len =n_len;
1066                 A->max_n_seq    =n_nseq;
1067
1068
1069                 A->seq_comment=new_realloc_char ( A->seq_comment, n_nseq, -1);
1070                 A->aln_comment=new_realloc_char ( A->aln_comment, n_nseq, -1);
1071
1072                 A->name   =new_realloc_char ( A->name, n_nseq, -1);
1073
1074
1075                 A->file   =new_realloc_char ( A->file, n_nseq, -1);
1076
1077                 A->tree_order   =new_realloc_char ( A->tree_order, n_nseq, -1);
1078                 A->seq_al =new_realloc_char ( A->seq_al, n_nseq, n_len);
1079                 A->order  =new_realloc_int  ( A->order, n_nseq, -1);
1080
1081                 if ( A->seq_cache) A->seq_cache=new_realloc_int  ( A->seq_cache, n_nseq,n_len);
1082                 if ( A->cdna_cache)A->cdna_cache=new_realloc_int  ( A->cdna_cache, n_nseq,n_len);
1083
1084
1085                 A->score_seq    =vrealloc( A->score_seq, sizeof (int)*(n_nseq));
1086                 for ( a=nseq; a< n_nseq; a++)A->score_seq[a]=0;
1087
1088
1089
1090             }
1091         return A;
1092         }
1093
1094
1095 long aln_stack (Alignment *A, int mode)
1096 {
1097   static long *list;
1098   static int size;
1099   static int max_size;
1100
1101
1102   if (A==NULL) return 0;
1103   else if ( mode==DECLARE_ALN)
1104     {
1105       if ( size==max_size)
1106         {
1107           max_size+=1000;
1108           list=vrealloc (list, max_size*sizeof (long));
1109         }
1110       list[size++]=(long)A;
1111       return 0;
1112     }
1113   else if (mode==FREE_ALN)
1114     {
1115       int a, b;
1116       for (a=0; a<size; a++)
1117         {
1118           if (list[a]==(long)A)
1119             {
1120               for (b=a+1; b<size;b++)
1121                 {
1122                   list[b-1]=list[b];
1123                 }
1124               list[b-1]=0;
1125               size--;
1126               return 1;
1127             }
1128         }
1129       return 0;
1130     }
1131   else if ( mode==EXTRACT_ALN)
1132     {
1133       return list[size--];
1134     }
1135   else
1136     {
1137       printf_exit (EXIT_FAILURE, stderr, "ERROR: Unknown mode for aln_stack");
1138       return 0;
1139     }
1140 }
1141
1142 Sequence* free_aln ( Alignment *LA){ return free_Alignment(LA);}
1143 Alignment *free_data_in_aln (Alignment *A)
1144 {
1145   //Frees only the sequence data (keeps profile information)
1146   A->seq_al=free_char (A->seq_al, -1);
1147   A->seq_comment=free_char (A->seq_comment, -1);
1148   A->aln_comment=free_char (A->aln_comment, -1);
1149   A->name=free_char (A->name, -1);
1150   A->expanded_order=free_char (A->expanded_order, -1);
1151
1152   A->order=free_int (A->order, -1);
1153   A->seq_cache=free_int (A->seq_cache, -1);
1154   A->cdna_cache=free_int (A->cdna_cache, -1);
1155   A->score_res=free_int (A->score_res, -1);
1156   free_sequence (A->S, -1);
1157   A->S=NULL;
1158   return A;
1159
1160 }
1161 Sequence* free_Alignment ( Alignment *LA)
1162         {
1163           /* Does not free the A->S field (sequences of A)*/
1164
1165           Sequence *S;
1166           //aln_stack checks the alignment has not already been freed
1167           if ( LA==NULL || !aln_stack(LA,FREE_ALN)){return NULL;}
1168
1169           S=LA->S;
1170           free_char ( LA->file, -1);
1171           free_char ( LA->seq_al, -1);
1172           free_int  ( LA->seq_cache, -1);
1173           free_int  ( LA->cdna_cache, -1);
1174           free_char ( LA->name,-1);
1175
1176           free_char ( LA->tree_order,-1);
1177           vfree ( LA->generic_comment);
1178           free_char ( LA->seq_comment, -1);
1179           free_char ( LA->aln_comment, -1);
1180
1181           free_int  ( LA->order, -1);
1182
1183           vfree ( LA->score_seq);
1184           vfree ( LA->len);
1185
1186           free_profile (LA->P);
1187           if ( LA->A){free_Alignment (LA->A);LA->A=NULL;}
1188
1189
1190           vfree ( LA);
1191           return S;
1192         }
1193
1194 Alignment * update_aln_random_tag ( Alignment *A)
1195 {
1196   static int tag;
1197   if ( !A) return A;
1198
1199   A->random_tag=++tag;
1200   return A;
1201 }
1202
1203 Profile   *copy_profile   (Profile *P1)
1204 {
1205
1206   Profile *P;
1207
1208   if ( !P1) return NULL;
1209   P=declare_profile ( P1->alphabet, P1->max_len);
1210   P->count=copy_int (P1->count, P->count, -1, -1);
1211   P->count2=copy_int (P1->count2, P->count2, -1, -1);
1212   P->count3=copy_int (P1->count3, P->count3, -1, -1);
1213
1214   return P;
1215
1216 }
1217
1218
1219 Profile   *declare_profile(char *alphabet, int len)
1220 {
1221   Profile *P;
1222   P=vcalloc ( 1, sizeof ( Profile));
1223   P->alp_size=strlen(alphabet);
1224   P->max_len=len;
1225   P->alphabet=vcalloc ( strlen (alphabet)+2, sizeof (char));
1226   sprintf ( P->alphabet, "%s", alphabet);
1227
1228   P->count=declare_int( P->alp_size+2, len);
1229   P->count2=declare_int(100, len);
1230   P->count3=declare_int(100, len);
1231
1232   return P;
1233 }
1234 Profile * free_profile ( Profile *P)
1235 {
1236   if (!P) return NULL;
1237   else
1238     {
1239       vfree (P->alphabet);
1240       free_int ( P->count, -1);
1241       free_int ( P->count2, -1);
1242       vfree (P);
1243     }
1244   return NULL;
1245 }
1246
1247
1248 /************************************************************************/
1249 /*                                                                      */
1250 /*             ALLOCATION                                               */
1251 /*                                                                      */
1252 /*                                                                      */
1253 /************************************************************************/
1254
1255
1256 double alloc_mem;
1257 double max_mem;
1258 double tot_mem;
1259 Memcontrol *memlast;
1260
1261 FILE* print_mem_usage (FILE *fp, char *comment)
1262 {
1263   fprintf ( fp, "# %s Memory Usage: Current= %.3f Mb, Max= %.3f Mb\n", comment,(float)((float)alloc_mem/(1024*1024)),(float)((float)tot_mem/(1024*1024)) );
1264   return fp;
1265 }
1266 void set_max_mem (int m)
1267 {
1268   max_mem=m*1024*1024;
1269 }
1270
1271 int verify_memory (int s)
1272 {
1273   alloc_mem+=s;
1274
1275   tot_mem=(alloc_mem>tot_mem)?alloc_mem:tot_mem;
1276
1277   if (max_mem && alloc_mem>max_mem)
1278     {
1279       fprintf (stderr, "\n%s Requires Too Much Memory: %d Megabytes [FATAL:%s]\n", PROGRAM,(int)(alloc_mem/1024*1024),PROGRAM);
1280       fprintf (stderr, "Tip: Rerun your Job with a smaller dataset\n");
1281
1282       myexit (EXIT_FAILURE);
1283     }
1284   else
1285     return 1;
1286   return 0;
1287
1288 }
1289
1290 int my_assert ( void *p, int index)
1291 {
1292   static int warning;
1293
1294   if (!warning)
1295     {
1296       fprintf ( stderr, "\n****************************************************************\n");
1297       fprintf ( stderr, "\n          DEBUG MODE [Rebuild For Better Performances]          \n");
1298       fprintf ( stderr, "\n*****************************************************************\n");
1299       warning=1;
1300     }
1301
1302   if ( !is_dynamic_memory(p)) return 1;
1303   else if ( read_array_size_new (p)<=index)
1304     {
1305       fprintf ( stderr, "\nFaulty Allocation: Size=%d Access=%d\n", read_array_size (p,0),index);
1306       return 0;
1307     }
1308   else
1309     {
1310       return 1;
1311     }
1312 }
1313
1314
1315
1316 void * vmalloc ( size_t size)
1317         {
1318         void * x;
1319         Memcontrol *M;
1320
1321         verify_memory (size+2*sizeof (Memcontrol));
1322
1323         if ( size==0)
1324             return NULL; /*crash ("\n0 bytes in vmalloc\n");*/
1325         else
1326             {
1327
1328               x= malloc (size + 2*sizeof (Memcontrol));
1329               //x=dlmalloc (size + 2*sizeof (Memcontrol));
1330             if ( x==NULL)
1331                 {
1332                   printf_exit (EXIT_FAILURE,stderr, "\nFAILED TO ALLOCATE REQUIRED MEMORY (vmalloc)\n");
1333
1334                 }
1335             else
1336               {
1337                 M=x;
1338                 M[0].size=size;
1339                 M[0].size_element=0;
1340                 sprintf ( M[0].check, "dy");
1341                 M+=2;
1342                 x=M;
1343                 return x;
1344               }
1345             }
1346         return NULL;}
1347
1348
1349
1350 void *vcalloc (size_t nobj, size_t size)
1351 {
1352   return sub_vcalloc (nobj,size, MEMSET0);
1353 }
1354 void *vcalloc_nomemset ( size_t nobj, size_t size)
1355 {
1356   return sub_vcalloc (nobj, size, NO_MEMSET0);
1357 }
1358 void *sub_vcalloc ( size_t nobj, size_t size, int MODE)
1359         {
1360         void *x;
1361         Memcontrol *M;
1362
1363         if ( nobj<=0 || size<=0)return NULL;/*crash ("\n0 bytes in vmalloc\n");*/
1364         else x=vmalloc (nobj*size);
1365
1366
1367         M=x;M-=2;M[0].size_element=size;M+=2;x=M;
1368
1369         if ( x==NULL)
1370                 {
1371                 crash ( "\nFAILED TO ALLOCATE REQUIRED MEMORY (vcalloc)\n");
1372                 return NULL;
1373                 }
1374         else
1375           {
1376             if ( MODE==MEMSET0)
1377               {
1378                 x=memset (x,0, nobj*size);
1379               }
1380             else
1381               {
1382                 if (nobj)x=memset (x, 0, size);
1383               }
1384             return x;
1385           }
1386         }
1387
1388 void *vrealloc ( void *p, size_t size)
1389         {
1390         void *x;
1391         Memcontrol *M;
1392         size_t i_size;
1393         int a;
1394
1395
1396         if ( p==NULL)
1397           {
1398             x=vmalloc (size);
1399             memset (x, 0, size);
1400
1401             return x;
1402           }
1403         else
1404           {
1405             M=p;
1406             M-=2;
1407             i_size=M[0].size;
1408             p=M;
1409
1410
1411             if ( size<=0){return NULL;vfree (p);return NULL;}
1412             else
1413               {
1414                 verify_memory (size - i_size);
1415                 x=realloc ( p, size+2*sizeof(Memcontrol));
1416
1417                 if ( x==NULL){crash ( "\nFAILED TO ALLOCATE REQUIRED MEMORY (realloc)\n");return NULL;}
1418                 M=x;
1419                 M[0].size=size;
1420                 M+=2;
1421                 x=M;
1422                 for ( a=i_size; a< size; a++)((char*)x)[a]=0;
1423                 return x;
1424               }
1425           }
1426         return NULL;
1427         }
1428 void vfree ( void *p)
1429      {
1430        Memcontrol *M;
1431        size_t size;
1432
1433        if ( !p)return;
1434        else
1435          {
1436            M=p;
1437            M-=2;
1438            size=M[0].size;
1439
1440            p=M;
1441            free(p);
1442
1443            verify_memory (-(size+2*sizeof(Memcontrol)));
1444          }
1445      }
1446 void vfree_all (void *p)
1447 {
1448   Memcontrol *n;
1449   while (memlast)
1450     {
1451       n=memlast->p;
1452       vfree (memlast+2);
1453       memlast=n;
1454     }
1455 }
1456 /*********************************************************************/
1457 /*                                                                   */
1458 /*                          SIZES                                    */
1459 /*                                                                   */
1460 /*                                                                   */
1461 /*********************************************************************/
1462 #define WRITE_SIZE(type,function)\
1463 void function ( int x, type *array, int os)\
1464      {\
1465      fprintf(stderr, "\nwrite_size is a deprecated function [Warning:%s]\n", PROGRAM);return;\
1466      }
1467 WRITE_SIZE(short,write_size_short)
1468 WRITE_SIZE(char,write_size_char)
1469 WRITE_SIZE(int,write_size_int)
1470 WRITE_SIZE(float,write_size_float)
1471 WRITE_SIZE(double,write_size_double)
1472
1473 #define READ_ARRAY_SIZE(type, function)\
1474 int function (void *array, size_t size)\
1475     {\
1476       return read_array_size (array, size);\
1477     }
1478 READ_ARRAY_SIZE(short,read_size_short)
1479 READ_ARRAY_SIZE(char,read_size_char)
1480 READ_ARRAY_SIZE(int,read_size_int)
1481 READ_ARRAY_SIZE(float,read_size_float)
1482 READ_ARRAY_SIZE(double,read_size_double)
1483
1484
1485 int read_array_size_new (void *array)
1486 {
1487   return read_array_size ( array, 0);
1488 }
1489 int read_array_size (void *array, size_t size)
1490     {
1491     Memcontrol *p;
1492     if (array==NULL)return 0;
1493     p=(Memcontrol *)array;
1494     p-=2;
1495     if ( p[0].size_element ==0 && size==0)
1496       {
1497         fprintf ( stderr, "\nERROR in read_array_size: trying to read the size of a malloced block");
1498       }
1499     else if ( size ==0) return (int)p[0].size/p[0].size_element;
1500
1501     return (int)p[0].size/size;
1502
1503     }
1504 int is_dynamic_memory ( void *array)
1505 {
1506   Memcontrol *p;
1507   if (array==NULL)return 0;
1508   p=(Memcontrol *)array;
1509   p-=2;
1510   if ( strm (p[0].check, "dy"))return 1;
1511   return 0;
1512 }
1513
1514 /************************************************************************/
1515 /*                                                                      */
1516 /*             DECLARE 2d ARRAYS                                        */
1517 /*                                                                      */
1518 /*                                                                      */
1519 /************************************************************************/
1520
1521 void * free_arrayN(void *p, int n)
1522 {
1523   int a, s;
1524   void **i;
1525
1526
1527   if ( p==NULL) return NULL;
1528   else if ( n==1)vfree ((void *)p);
1529   else
1530     {
1531       i=(void**)p;
1532       s=read_array_size ( (void *)p, sizeof ( void *));
1533       for ( a=0; a< s; a++)free_arrayN ((void *)i[a], n-1);
1534       vfree (p);
1535     }
1536   return NULL;
1537 }
1538
1539 void * declare_arrayNnomemset (int ndim, size_t size, ...)
1540 {
1541   va_list ap;
1542   int *array;
1543   void **p;
1544   int a;
1545
1546   va_start (ap, size);
1547
1548   array=vcalloc (ndim, sizeof (int));
1549   for ( a=0; a< ndim; a++)
1550     {
1551       array[a]=va_arg (ap,int);
1552       if ( array[a]<0){va_end(ap);return NULL;}
1553
1554     }
1555   va_end (ap);
1556
1557   if ( ndim==2)
1558     {
1559
1560       p=vcalloc_nomemset (array[0], sizeof ( void*));
1561       for (a=0; a< array[0]; a++)
1562         {
1563         p[a]=vcalloc_nomemset (array[1], size);
1564         }
1565     }
1566   else
1567     {
1568       p=declare_arrayN2nomemset (ndim, array, size);
1569     }
1570   vfree (array);
1571   return p;
1572 }
1573
1574 void *declare_arrayN2nomemset ( int ndim, int *A, size_t size)
1575 {
1576   int a;
1577   void **p;
1578
1579   if ( ndim>1)
1580     {
1581       p=vcalloc_nomemset (A[0], sizeof (void*));
1582       for ( a=0; a<A[0]; a++)
1583         p[a]=declare_arrayN2(ndim-1, A+1, size);
1584     }
1585   else
1586     {
1587       p=vcalloc_nomemset (A[0], size);
1588     }
1589   return p;
1590 }
1591
1592 void * declare_arrayN (int ndim, size_t size, ...)
1593 {
1594   va_list ap;
1595   int *array;
1596   void **p;
1597   int a;
1598
1599
1600   va_start (ap, size);
1601
1602   array=vcalloc (ndim, sizeof (int));
1603   for ( a=0; a< ndim; a++)
1604     {
1605       array[a]=va_arg (ap,int);
1606       if ( array[a]<0){va_end(ap);return NULL;}
1607
1608     }
1609   va_end (ap);
1610
1611   if ( ndim==2)
1612     {
1613
1614       p=vcalloc_nomemset (array[0], sizeof ( void*));
1615       for (a=0; a< array[0]; a++)
1616         p[a]=vcalloc (array[1], size);
1617     }
1618   else
1619     {
1620       p=declare_arrayN2 (ndim, array, size);
1621     }
1622   vfree (array);
1623   return p;
1624 }
1625
1626 void *declare_arrayN2 ( int ndim, int *A, size_t size)
1627 {
1628   int a;
1629   void **p;
1630
1631   if ( ndim>1)
1632     {
1633       p=vcalloc_nomemset (A[0], sizeof (void*));
1634       for ( a=0; a<A[0]; a++)
1635         p[a]=declare_arrayN2(ndim-1, A+1, size);
1636     }
1637   else
1638     {
1639       p=vcalloc (A[0], size);
1640     }
1641   return p;
1642 }
1643
1644 void **declare_array (int first, int second, size_t size)
1645 {
1646   return (void **)declare_arrayN (2,size,first, second);
1647 }
1648
1649
1650 #define DECLARE_ARRAY(type,wf,rf,function)\
1651 type**  function (int first, int second)\
1652   {\
1653     return (type **)declare_arrayN (2,sizeof(type), first, second);\
1654    }
1655
1656 int **declare_int2 (int f, int *s, int d)
1657 {
1658   int **r;
1659   int a;
1660   r=vcalloc ( f, sizeof (int*));
1661   for (a=0; a<f; a++)
1662     r[a]=vcalloc (s[a]+d, sizeof (int));
1663   return r;
1664 }
1665
1666
1667
1668
1669
1670 DECLARE_ARRAY(short,write_size_short,read_size_short,declare_short)
1671 DECLARE_ARRAY(char,write_size_char,read_size_char,declare_char)
1672 DECLARE_ARRAY(int,write_size_int,read_size_int,declare_int)
1673 DECLARE_ARRAY(float,write_size_float,read_size_float,declare_float)
1674 DECLARE_ARRAY(double,write_size_double,read_size_double,declare_double)
1675
1676 void **declare_array_nomemset (int first, int second, size_t size)
1677 {
1678   return (void **)declare_arrayNnomemset (2,size,first, second);
1679 }
1680 #define DECLARE_ARRAY_NMS(type,wf,rf,function)\
1681 type**  function (int first, int second)\
1682   {\
1683     return (type **)declare_arrayNnomemset (2,sizeof(type), first, second);\
1684    }
1685 DECLARE_ARRAY_NMS(short,write_size_short,read_size_short,declare_short_nomemset)
1686 DECLARE_ARRAY_NMS(char,write_size_char,read_size_char,declare_char_nomemset)
1687 DECLARE_ARRAY_NMS(int,write_size_int,read_size_int,declare_int_nomemset)
1688 DECLARE_ARRAY_NMS(float,write_size_float,read_size_float,declare_float_nomemset)
1689 DECLARE_ARRAY_NMS(double,write_size_double,read_size_double,declare_double_nomemset)
1690
1691
1692 Alignment ** declare_aln_array ( int first)
1693     {
1694     Alignment ** array;
1695     int a;
1696
1697
1698     array=vcalloc (first, sizeof (Alignment*));
1699     for ( a=0; a< first; a++)
1700       array[a]=declare_Alignment (NULL);
1701     return array;
1702     }
1703
1704
1705
1706
1707 /************************************************************************/
1708 /*                                                                      */
1709 /*             Realloc 2d ARRAYS                                        */
1710 /*                                                                      */
1711 /*                                                                      */
1712 /************************************************************************/
1713 void ** realloc_arrayN(int ndim,void **main_array,size_t size, ...)
1714 {
1715   va_list ap;
1716   int *array;
1717   void **p;
1718   int a;
1719
1720   /*dim size==-1: keep current size*/
1721   /*Dim sizes are the absolute size (not the extension*/
1722   /*If array getting shorter, memory is Not Claimed back*/
1723
1724   array=vcalloc (ndim, sizeof (int));
1725   va_start (ap, size);
1726   for ( a=0; a< ndim; a++)
1727     {
1728       array[a]=va_arg (ap,int);
1729       if ( array[a]<-1){va_end(ap);return NULL;}
1730     }
1731   va_end (ap);
1732
1733   p=realloc_arrayN2 (ndim, main_array, array,size);
1734   vfree (array);
1735   return p;
1736 }
1737
1738 void **realloc_arrayN2 ( int ndim, void ** p, int *A, size_t size)
1739 {
1740   int a;
1741   int o, n;
1742
1743
1744   if (ndim==0)return NULL;
1745
1746   if ( ndim>1)
1747     {
1748       o=read_array_size (p,sizeof (void*));
1749       if (A[0]>o)p=vrealloc (p, sizeof (void*)*A[0]);
1750       n=(A[0]==-1)?o:A[0];
1751       for ( a=0; a<n; a++)
1752         p[a]=realloc_arrayN2(ndim-1,p[a], A+1, size);
1753     }
1754   else
1755     {
1756       o=read_array_size (p, size);
1757       if (A[0]>o)p=vrealloc (p, size*A[0]);
1758     }
1759   return p;
1760 }
1761
1762
1763
1764 void ** realloc_array (void **array,size_t size, int first, int second, int ext1, int ext2)
1765 {
1766   int a;
1767   int d1, d2;
1768   if ( array==NULL)return declare_array (((first==-1)?0:first)+ext1, ((second==-1)?0:second)+ext2, size);
1769   else if ( first==-1)
1770     {
1771       first=read_array_size (array, sizeof (void *));
1772     }
1773   if (second==-1)second=read_array_size(array[0], size);
1774
1775   d1=first+ext1;
1776   d2=second+ext2;
1777
1778   for ( a=d1; a<first; a++)vfree (array[a]);
1779   array=vrealloc (array, sizeof (void*)*d1);
1780
1781   if ( d2!=second)
1782     {
1783       for (a=0; a<d1 && a<first; a++)
1784         array[a]=vrealloc ( array[a], size*d2);
1785     }
1786   for ( a=first; a< d1; a++)
1787     array[a]=vrealloc ( array[a], size*d2);
1788   return array;
1789 }
1790
1791
1792
1793
1794
1795 #define REALLOC_ARRAY(type,wf,rf,function1,function2,function3)\
1796 type ** function1 ( type **array, int first, int second, int ext1, int ext2)\
1797     {\
1798 return (type **)realloc_array ((void **)array, sizeof (type),first, second, ext1, ext2);\
1799 \
1800     }
1801 REALLOC_ARRAY(short,write_size_short,read_size_short,realloc_short,declare_short,free_short)
1802 REALLOC_ARRAY(char,write_size_char,read_size_char,realloc_char,declare_char,free_char)
1803 REALLOC_ARRAY(int,write_size_int,read_size_int,realloc_int,declare_int,free_int)
1804 REALLOC_ARRAY(float,write_size_float,read_size_float,realloc_float,declare_float,free_float)
1805 REALLOC_ARRAY(double,write_size_double,read_size_double,realloc_double,declare_double,free_double)
1806
1807
1808
1809 #define NEW_REALLOC_ARRAY(type,wf,rf,function1,function2,function3)\
1810 type ** function1 ( type **array, int ext1, int ext2)\
1811     {int a, b;\
1812      int first, l1;\
1813      int second, l2;\
1814      type **new_array;\
1815 \
1816      first=rf(array,sizeof (type*));\
1817      second=rf(array[0],sizeof (type));\
1818      \
1819      if ( ext1==-1)ext1=first;\
1820      if ( ext2==-1)ext2=second;\
1821      l1=MIN(ext1, first);\
1822      l2=MIN(ext2, second);\
1823      new_array=declare_arrayN(2,sizeof(type),ext1, ext2);\
1824      for ( a=0; a<l1; a++)\
1825        for ( b=0; b<l2; b++)new_array[a][b]=array[a][b];\
1826      function3(array, -1);\
1827      return new_array;\
1828     }
1829
1830 NEW_REALLOC_ARRAY(short,write_size_short,read_size_short,new_realloc_short,declare_short,free_short)
1831 NEW_REALLOC_ARRAY(char,write_size_char,read_size_char,new_realloc_char,declare_char,free_char)
1832 NEW_REALLOC_ARRAY(int,write_size_int,read_size_int,new_realloc_int,declare_int,free_int)
1833 NEW_REALLOC_ARRAY(float,write_size_float,read_size_float,new_realloc_float,declare_float,free_float)
1834 NEW_REALLOC_ARRAY(double,write_size_double,read_size_double,new_realloc_double,declare_double,free_double)
1835 Alignment ** realloc_aln_array ( Alignment **array, int ext1)
1836     {
1837     int a;
1838     int first;
1839
1840
1841     if ( array==NULL)
1842          {
1843          array=declare_aln_array(ext1);
1844          return array;
1845          }
1846      first=read_array_size ( array, sizeof ( Alignment *));
1847      if ( ext1>0)
1848         {
1849         array=vrealloc ( array, (sizeof (Alignment*))*(first+ext1));
1850         for ( a=first; a<first+ext1; a++)array[a]=declare_Alignment (NULL);
1851         }
1852     else if ( ext1==0);
1853     else if ( ext1<0)
1854          {
1855          for ( a=first-1; a>=(first+ext1);a--)free_Alignment (array[a]);
1856          array=vrealloc ( array, (sizeof (Alignment*))*(first+ext1));
1857          }
1858      return array;
1859     }
1860
1861 /************************************************************************/
1862 /*                                                                      */
1863 /*            free 2d ARRAYS                                        */
1864 /*                                                                      */
1865 /*                                                                      */
1866 /************************************************************************/
1867 #define FREE_ARRAY(type,wf,rf,function) \
1868 type ** function (type **array, int first)\
1869     {\
1870       return free_arrayN((void*)array, 2);\
1871     }
1872 FREE_ARRAY(short,write_size_short,read_size_short,free_short)
1873 FREE_ARRAY(char,write_size_char,read_size_char,free_char)
1874 FREE_ARRAY(int,write_size_int,read_size_int,free_int)
1875 FREE_ARRAY(float,write_size_float,read_size_float,free_float)
1876 FREE_ARRAY(double,write_size_double,read_size_double,free_double)
1877
1878
1879
1880 Alignment ** free_aln_array (Alignment **array)
1881    {
1882    int a;
1883    int len;
1884
1885
1886    if ( array==NULL)return NULL;
1887    len=read_array_size ( array, sizeof (Alignment *));
1888    for ( a=1; a< len; a++)free_Alignment(array[a]);
1889    vfree ( array);
1890    return  NULL;
1891    }
1892
1893 Fname *declare_fname (int size)
1894    {
1895    Fname *F;
1896
1897    size+=strlen (get_home_4_tcoffee())+FILENAMELEN+1;
1898
1899    F=vcalloc ( 1, sizeof (Fname));
1900    F->name  =vcalloc ( size, sizeof (char));
1901    F->path  =vcalloc ( size, sizeof (char));
1902    F->suffix=vcalloc ( size, sizeof (char));
1903    F->full=vcalloc ( size, sizeof (char));
1904    return F;
1905    }
1906
1907 Fname *free_fname ( Fname *F)
1908    {
1909    vfree (F->name);
1910    vfree (F->path);
1911    vfree (F->suffix);
1912    return NULL;
1913    }
1914 /******************************COPYRIGHT NOTICE*******************************/
1915 /*© Centro de Regulacio Genomica */
1916 /*and */
1917 /*Cedric Notredame */
1918 /*Fri Feb 18 08:27:45 CET 2011 - Revision 596. */
1919 /*All rights reserved.*/
1920 /*This file is part of T-COFFEE.*/
1921 /**/
1922 /*    T-COFFEE is free software; you can redistribute it and/or modify*/
1923 /*    it under the terms of the GNU General Public License as published by*/
1924 /*    the Free Software Foundation; either version 2 of the License, or*/
1925 /*    (at your option) any later version.*/
1926 /**/
1927 /*    T-COFFEE is distributed in the hope that it will be useful,*/
1928 /*    but WITHOUT ANY WARRANTY; without even the implied warranty of*/
1929 /*    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the*/
1930 /*    GNU General Public License for more details.*/
1931 /**/
1932 /*    You should have received a copy of the GNU General Public License*/
1933 /*    along with Foobar; if not, write to the Free Software*/
1934 /*    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA*/
1935 /*...............................................                                                                                      |*/
1936 /*  If you need some more information*/
1937 /*  cedric.notredame@europe.com*/
1938 /*...............................................                                                                                                                                     |*/
1939 /**/
1940 /**/
1941 /*      */
1942 /******************************COPYRIGHT NOTICE*******************************/