JAL-1681 show cDNA consensus on protein alignment - first version
[jalview.git] / src / jalview / analysis / StructureFrequency.java
index cae5da2..cb78b04 100644 (file)
@@ -1,28 +1,33 @@
 /*
- * Jalview - A Sequence Alignment Editor and Viewer (Version 2.8)
- * Copyright (C) 2012 J Procter, AM Waterhouse, LM Lui, J Engelhardt, G Barton, M Clamp, S Searle
+ * Jalview - A Sequence Alignment Editor and Viewer (Version 2.8.2)
+ * Copyright (C) 2014 The Jalview Authors
  * 
  * This file is part of Jalview.
  * 
  * Jalview is free software: you can redistribute it and/or
  * modify it under the terms of the GNU General Public License 
- * as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
+ * as published by the Free Software Foundation, either version 3
+ * of the License, or (at your option) any later version.
  *  
  * Jalview is distributed in the hope that it will be useful, but 
  * WITHOUT ANY WARRANTY; without even the implied warranty 
  * of MERCHANTABILITY or FITNESS FOR A PARTICULAR 
  * PURPOSE.  See the GNU General Public License for more details.
  * 
- * You should have received a copy of the GNU General Public License along with Jalview.  If not, see <http://www.gnu.org/licenses/>.
+ * You should have received a copy of the GNU General Public License
+ * along with Jalview.  If not, see <http://www.gnu.org/licenses/>.
+ * The Jalview Authors are detailed in the 'AUTHORS' file.
  */
-
-
 package jalview.analysis;
 
-import java.util.*;
-
+import jalview.datamodel.AlignmentAnnotation;
+import jalview.datamodel.Annotation;
+import jalview.datamodel.SequenceFeature;
+import jalview.datamodel.SequenceI;
 import jalview.util.Format;
-import jalview.datamodel.*;
+
+import java.util.ArrayList;
+import java.util.Hashtable;
 
 /**
  * Takes in a vector or array of sequences and column start and column end and
@@ -35,6 +40,8 @@ import jalview.datamodel.*;
  */
 public class StructureFrequency
 {
+  public static final int STRUCTURE_PROFILE_LENGTH = 74;
+
   // No need to store 1000s of strings which are not
   // visible to the user.
   public static final String MAXCOUNT = "C";
@@ -60,15 +67,15 @@ public class StructureFrequency
    */
   public static int findPair(SequenceFeature[] pairs, int indice)
   {
-         System.out.print("indice"+indice+"    ");
+
     for (int i = 0; i < pairs.length; i++)
     {
       if (pairs[i].getBegin() == indice)
-        
+
       {
-         System.out.println(pairs[i].getEnd());
+
         return pairs[i].getEnd();
-        
+
       }
     }
     return -1;
@@ -89,19 +96,11 @@ public class StructureFrequency
           int end, Hashtable[] result, boolean profile,
           AlignmentAnnotation rnaStruc)
   {
-//     System.out.println("longueur="+sequences.length);
-//     for(int l=0;l<=(sequences.length-1);l++){  
-//     System.out.println("sequences "+l+":"+sequences[l].getSequenceAsString());
-//     }
-//     System.out.println("start="+start);
-       System.out.println("end="+end);
-//     System.out.println("result="+result.length);
-//
-//     System.out.println("profile="+profile);
-//     System.out.println("rnaStruc="+rnaStruc);
+
     Hashtable residueHash;
     String maxResidue;
     char[] struc = rnaStruc.getRNAStruc().toCharArray();
+
     SequenceFeature[] rna = rnaStruc._rnasecstr;
     char c, s, cEnd;
     int count = 0, nonGap = 0, i, bpEnd = -1, j, jSize = sequences.length;
@@ -109,7 +108,6 @@ public class StructureFrequency
     int[][] pairs;
     float percentage;
     boolean wooble = true;
-
     for (i = start; i < end; i++) // foreach column
     {
       residueHash = new Hashtable();
@@ -117,11 +115,11 @@ public class StructureFrequency
       values = new int[255];
       pairs = new int[255][255];
       bpEnd = -1;
-      //System.out.println("s="+struc[i]);
+      // System.out.println("s="+struc[i]);
       if (i < struc.length)
       {
         s = struc[i];
-        
+
       }
       else
       {
@@ -141,76 +139,74 @@ public class StructureFrequency
       }
       else
       {
-        
-         
+
         bpEnd = findPair(rna, i);
-       
-        if (bpEnd>-1)
-        {
-        for (j = 0; j < jSize; j++) // foreach row
+
+        if (bpEnd > -1)
         {
-          if (sequences[j] == null)
+          for (j = 0; j < jSize; j++) // foreach row
           {
-            System.err
-                    .println("WARNING: Consensus skipping null sequence - possible race condition.");
-            continue;
-          }
-          c = sequences[j].getCharAt(i);
-          //System.out.println("c="+c);
-          
+            if (sequences[j] == null)
+            {
+              System.err
+                      .println("WARNING: Consensus skipping null sequence - possible race condition.");
+              continue;
+            }
+            c = sequences[j].getCharAt(i);
+            // System.out.println("c="+c);
 
             // standard representation for gaps in sequence and structure
             if (c == '.' || c == ' ')
             {
-              System.err
-                      .println("WARNING: Consensus skipping null sequence - possible race condition.");
+              c = '-';
+            }
+
+            if (c == '-')
+            {
+              values['-']++;
               continue;
             }
             cEnd = sequences[j].getCharAt(bpEnd);
-         
-         
-            System.out.println("pairs ="+c+","+cEnd);
-            if (checkBpType(c, cEnd)==true)
+
+            // System.out.println("pairs ="+c+","+cEnd);
+            if (checkBpType(c, cEnd) == true)
             {
               values['(']++; // H means it's a helix (structured)
               maxResidue = "(";
-              wooble=true;
-              System.out.println("It's a pair wc");
-              
+              wooble = true;
+              // System.out.println("It's a pair wc");
+
             }
-            if (checkBpType(c, cEnd)==false)
+            if (checkBpType(c, cEnd) == false)
             {
-              wooble =false;
+              wooble = false;
               values['[']++; // H means it's a helix (structured)
               maxResidue = "[";
-              System.out.println("It's an pair non canonic");
-              System.out.println(sequences[j].getRNA());
-              System.out.println(rnaStruc.getRNAStruc().charAt(i));
+
             }
-            pairs[c][cEnd]++;  
-           
-          
-               }
+            pairs[c][cEnd]++;
 
+          }
         }
         // nonGap++;
       }
       // UPDATE this for new values
       if (profile)
       {
+        // TODO 1-dim array with jsize in [0], nongapped in [1]; or Pojo
         residueHash.put(PROFILE, new int[][]
         { values, new int[]
         { jSize, (jSize - values['-']) } });
 
         residueHash.put(PAIRPROFILE, pairs);
       }
-      if (wooble==true)
+      if (wooble == true)
       {
-      count = values['('];
+        count = values['('];
       }
-      if (wooble==false)
+      if (wooble == false)
       {
-      count = values['['];
+        count = values['['];
       }
       residueHash.put(MAXCOUNT, new Integer(count));
       residueHash.put(MAXRESIDUE, maxResidue);
@@ -226,18 +222,20 @@ public class StructureFrequency
       }
       if (bpEnd > 0)
       {
-       values[')'] = values['('];
+        values[')'] = values['('];
         values[']'] = values['['];
         values['('] = 0;
         values['['] = 0;
         residueHash = new Hashtable();
-        if (wooble==true){
-               System.out.println(maxResidue+","+wooble);
-               maxResidue = ")";
+        if (wooble == true)
+        {
+          // System.out.println(maxResidue+","+wooble);
+          maxResidue = ")";
         }
-        if(wooble==false){
-               System.out.println(maxResidue+","+wooble);
-               maxResidue = "]";
+        if (wooble == false)
+        {
+          // System.out.println(maxResidue+","+wooble);
+          maxResidue = "]";
         }
         if (profile)
         {
@@ -255,6 +253,7 @@ public class StructureFrequency
         residueHash.put(PID_GAPS, new Float(percentage));
 
         result[bpEnd] = residueHash;
+
       }
     }
   }
@@ -355,18 +354,19 @@ public class StructureFrequency
       // initialised properly
       return;
     }
-    String fmtstr="%3.1f";
-    int precision=2;
-    while (nseq>100) {
+    String fmtstr = "%3.1f";
+    int precision = 2;
+    while (nseq > 100)
+    {
       precision++;
-      nseq/=10;
+      nseq /= 10;
     }
-    if (precision>2)
+    if (precision > 2)
     {
-      fmtstr = "%"+(2+precision)+"."+precision+"f";
+      fmtstr = "%" + (2 + precision) + "." + precision + "f";
     }
     Format fmt = new Format(fmtstr);
-    
+
     for (int i = iStart; i < width; i++)
     {
       Hashtable hci;
@@ -443,7 +443,8 @@ public class StructureFrequency
             tval = (vl[c] * 100f / profile[1][ignoreGapsInConsensusCalculation ? 1
                     : 0]);
             mouseOver += ((p == 0) ? "" : "; ") + (char) ((int[]) ca[c])[0]
-                    + (char) ((int[]) ca[c])[1] + " " + fmt.form(tval) + "%";
+                    + (char) ((int[]) ca[c])[1] + " " + fmt.form(tval)
+                    + "%";
             p++;
 
           }
@@ -469,16 +470,19 @@ public class StructureFrequency
   public static int[] extractProfile(Hashtable hconsensus,
           boolean ignoreGapsInConsensusCalculation)
   {
-    int[] rtnval = new int[74]; // 2*(5*5)+2
+    int[] rtnval = new int[STRUCTURE_PROFILE_LENGTH]; // 2*(5*5)+2
     int[][] profile = (int[][]) hconsensus.get(StructureFrequency.PROFILE);
     int[][] pairs = (int[][]) hconsensus
             .get(StructureFrequency.PAIRPROFILE);
 
     if (profile == null)
+    {
       return null;
+    }
 
     // TODO fix the object length, also do it in completeConsensus
-    Object[] ca = new Object[625];
+    // Object[] ca = new Object[625];
+    int[][] ca = new int[625][];
     float[] vl = new float[625];
     int x = 0;
     for (int c = 65; c < 90; c++)
@@ -493,21 +497,28 @@ public class StructureFrequency
     }
     jalview.util.QuickSort.sort(vl, ca);
 
-    rtnval[0] = 2;
+    int valuesCount = 0;
     rtnval[1] = 0;
+    int offset = 2;
     for (int c = 624; c > 0; c--)
     {
       if (vl[c] > 0)
       {
-        rtnval[rtnval[0]++] = ((int[]) ca[c])[0];
-        rtnval[rtnval[0]++] = ((int[]) ca[c])[1];
-        rtnval[rtnval[0]] = (int) (vl[c] * 100f / profile[1][ignoreGapsInConsensusCalculation ? 1
+        rtnval[offset++] = ca[c][0];
+        rtnval[offset++] = ca[c][1];
+        rtnval[offset] = (int) (vl[c] * 100f / profile[1][ignoreGapsInConsensusCalculation ? 1
                 : 0]);
-        rtnval[1] += rtnval[rtnval[0]++];
+        rtnval[1] += rtnval[offset++];
+        valuesCount++;
       }
     }
+    rtnval[0] = valuesCount;
 
-    return rtnval;
+    // insert profile type code in position 0
+    int[] result = new int[rtnval.length + 1];
+    result[0] = AlignmentAnnotation.STRUCTURE_PROFILE;
+    System.arraycopy(rtnval, 0, result, 1, rtnval.length);
+    return result;
   }
 
   public static void main(String args[])