[Seqinr-commits] r1603 - pkg/src

[noreply at r-forge.r-project.org]

Author: lobry
Date: 2009-04-26 13:49:05 +0200 (Sun, 26 Apr 2009)
New Revision: 1603

Modified:
   pkg/src/kaks.c
Log:
now any non ACGT base in a codon force it to a gap codon ---

Modified: pkg/src/kaks.c
===================================================================
--- pkg/src/kaks.c	2009-04-26 11:46:35 UTC (rev 1602)
+++ pkg/src/kaks.c	2009-04-26 11:49:05 UTC (rev 1603)
@@ -152,33 +152,34 @@
 
 /******************************************************************************/
 /*                                                                            */
-/* Replace codons with ambiguous bases with ---                               */
+/* Replace codons with non ACGT bases with ---                                */
 /*                                                                            */
 /******************************************************************************/
 
-	for (i=0;i<totseqs;i++){
-		for(j=0;j<lgseq;j++){
-			if ((*(seqIn[i]+j)!='A') && (*(seqIn[i]+j)!='G') && (*(seqIn[i]+j)!='C') && (*(seqIn[i]+j)!='T') && (*(seqIn[i]+j)!='-') ) {
-				if (j%3==0) {
-					*(seqIn[i]+j)='-';
-					*(seqIn[i]+j+1)='-';
-					*(seqIn[i]+j+2)='-';
-				}
-
-				if (j%3==1) {
-					*(seqIn[i]+j)='-';
-					*(seqIn[i]+j+1)='-';
-					*(seqIn[i]+j-1)='-';
-				}
-
-				if (j%3==2) {
-					*(seqIn[i]+j)='-';
-					*(seqIn[i]+j-1)='-';
-					*(seqIn[i]+j-2)='-';
-				}
-			}
-		}
+  for (i = 0 ; i < totseqs ; i++){
+    for(j = 0 ; j < lgseq ; j++){
+      if ((*(seqIn[i] + j) != 'A') && (*(seqIn[i] + j) != 'G') && (*(seqIn[i] + j) != 'C') && (*(seqIn[i] + j) != 'T') ) {
+        /* Base in first codon position */
+        if (j % 3 == 0) {
+	  *(seqIn[i] + j) = '-';
+	  *(seqIn[i] + j + 1) = '-';
+	  *(seqIn[i] + j + 2) = '-';
 	}
+        /* Base in second codon position */
+	if (j % 3 == 1) {
+	  *(seqIn[i] + j) = '-';
+	  *(seqIn[i] + j + 1) = '-';
+	  *(seqIn[i] + j - 1) = '-';
+	}
+        /* Base in third codon position */
+	if (j % 3 == 2) {
+	  *(seqIn[i] + j) = '-';
+	  *(seqIn[i] + j - 1) = '-';
+	  *(seqIn[i] + j - 2) = '-';
+	}
+      }
+    }
+  }
 
 /******************************************************************************/
 /*                                                                            */
@@ -314,11 +315,11 @@
 
   sat = sat1 = sat2 = 2;
   /*
-     Internal check at C level: 
+     Internal check at C level: this should be no more be necessary, I'll keep it just in case. JRL - 26-APR-2009
   */
   flgseq = (double) lgseq;
   if (flgseq / trois != lgseq / 3) {
-    REprintf("Error: the number of nucleotide after gap removal is not a multiple of 3.\nDid you align the CDS at the aa level?\nSee reverse.align().\n");
+    REprintf("Fatal error: the number of nucleotide after gap removal is not a multiple of 3.\nPlease report this bug on the seqinr diffusion list.\n");
     return(0); /* Should be R's NA but an int is returned by fastlwl */
   }
   for (i = 0; i < nbseq - 1; i++) {