Clover icon

jalviewX

  1. Project Clover database Wed Oct 31 2018 15:13:58 GMT
  2. Package jalview.analysis

File AlignSeq.java

 

Coverage histogram

../../img/srcFileCovDistChart8.png
19% of files have more coverage

Code metrics

162
349
30
1
1,125
712
122
0.35
11.63
30
4.07

Classes

Class Line # Actions
AlignSeq 51 349 122 144
0.733826273.4%
 

Contributing tests

This file is covered by 208 tests. .

Source view

1    /*
2    * Jalview - A Sequence Alignment Editor and Viewer ($$Version-Rel$$)
3    * Copyright (C) $$Year-Rel$$ The Jalview Authors
4    *
5    * This file is part of Jalview.
6    *
7    * Jalview is free software: you can redistribute it and/or
8    * modify it under the terms of the GNU General Public License
9    * as published by the Free Software Foundation, either version 3
10    * of the License, or (at your option) any later version.
11    *
12    * Jalview is distributed in the hope that it will be useful, but
13    * WITHOUT ANY WARRANTY; without even the implied warranty
14    * of MERCHANTABILITY or FITNESS FOR A PARTICULAR
15    * PURPOSE. See the GNU General Public License for more details.
16    *
17    * You should have received a copy of the GNU General Public License
18    * along with Jalview. If not, see <http://www.gnu.org/licenses/>.
19    * The Jalview Authors are detailed in the 'AUTHORS' file.
20    */
21    package jalview.analysis;
22   
23    import jalview.analysis.scoremodels.PIDModel;
24    import jalview.analysis.scoremodels.ScoreMatrix;
25    import jalview.analysis.scoremodels.ScoreModels;
26    import jalview.analysis.scoremodels.SimilarityParams;
27    import jalview.datamodel.AlignmentAnnotation;
28    import jalview.datamodel.AlignmentI;
29    import jalview.datamodel.Mapping;
30    import jalview.datamodel.Sequence;
31    import jalview.datamodel.SequenceI;
32    import jalview.util.Comparison;
33    import jalview.util.Format;
34    import jalview.util.MapList;
35    import jalview.util.MessageManager;
36   
37    import java.awt.Color;
38    import java.awt.Graphics;
39    import java.io.PrintStream;
40    import java.util.ArrayList;
41    import java.util.Arrays;
42    import java.util.List;
43    import java.util.StringTokenizer;
44   
45    /**
46    *
47    *
48    * @author $author$
49    * @version $Revision$
50    */
 
51    public class AlignSeq
52    {
53    private static final int MAX_NAME_LENGTH = 30;
54   
55    private static final int GAP_OPEN_COST = 120;
56   
57    private static final int GAP_EXTEND_COST = 20;
58   
59    private static final int GAP_INDEX = -1;
60   
61    public static final String PEP = "pep";
62   
63    public static final String DNA = "dna";
64   
65    private static final String NEWLINE = System.lineSeparator();
66   
67    float[][] score;
68   
69    float[][] E;
70   
71    float[][] F;
72   
73    int[][] traceback; // todo is this actually used?
74   
75    int[] seq1;
76   
77    int[] seq2;
78   
79    SequenceI s1;
80   
81    SequenceI s2;
82   
83    public String s1str;
84   
85    public String s2str;
86   
87    int maxi;
88   
89    int maxj;
90   
91    int[] aseq1;
92   
93    int[] aseq2;
94   
95    public String astr1 = "";
96   
97    public String astr2 = "";
98   
99    /** DOCUMENT ME!! */
100    public int seq1start;
101   
102    /** DOCUMENT ME!! */
103    public int seq1end;
104   
105    /** DOCUMENT ME!! */
106    public int seq2start;
107   
108    public int seq2end;
109   
110    int count;
111   
112    public float maxscore;
113   
114    int prev = 0;
115   
116    StringBuffer output = new StringBuffer();
117   
118    String type; // AlignSeq.PEP or AlignSeq.DNA
119   
120    private ScoreMatrix scoreMatrix;
121   
122    /**
123    * Creates a new AlignSeq object.
124    *
125    * @param s1
126    * first sequence for alignment
127    * @param s2
128    * second sequence for alignment
129    * @param type
130    * molecule type, either AlignSeq.PEP or AlignSeq.DNA
131    */
 
132  179 toggle public AlignSeq(SequenceI s1, SequenceI s2, String type)
133    {
134  179 seqInit(s1, s1.getSequenceAsString(), s2, s2.getSequenceAsString(),
135    type);
136    }
137   
138    /**
139    * Creates a new AlignSeq object.
140    *
141    * @param s1
142    * DOCUMENT ME!
143    * @param s2
144    * DOCUMENT ME!
145    * @param type
146    * DOCUMENT ME!
147    */
 
148  2 toggle public AlignSeq(SequenceI s1, String string1, SequenceI s2,
149    String string2, String type)
150    {
151  2 seqInit(s1, string1.toUpperCase(), s2, string2.toUpperCase(), type);
152    }
153   
154    /**
155    * DOCUMENT ME!
156    *
157    * @return DOCUMENT ME!
158    */
 
159  22 toggle public float getMaxScore()
160    {
161  22 return maxscore;
162    }
163   
164    /**
165    * DOCUMENT ME!
166    *
167    * @return DOCUMENT ME!
168    */
 
169  460 toggle public int getSeq2Start()
170    {
171  460 return seq2start;
172    }
173   
174    /**
175    * DOCUMENT ME!
176    *
177    * @return DOCUMENT ME!
178    */
 
179  145 toggle public int getSeq2End()
180    {
181  145 return seq2end;
182    }
183   
184    /**
185    * DOCUMENT ME!
186    *
187    * @return DOCUMENT ME!
188    */
 
189  460 toggle public int getSeq1Start()
190    {
191  460 return seq1start;
192    }
193   
194    /**
195    * DOCUMENT ME!
196    *
197    * @return DOCUMENT ME!
198    */
 
199  145 toggle public int getSeq1End()
200    {
201  145 return seq1end;
202    }
203   
204    /**
205    * DOCUMENT ME!
206    *
207    * @return DOCUMENT ME!
208    */
 
209  2 toggle public String getOutput()
210    {
211  2 return output.toString();
212    }
213   
214    /**
215    * DOCUMENT ME!
216    *
217    * @return DOCUMENT ME!
218    */
 
219  148 toggle public String getAStr1()
220    {
221  148 return astr1;
222    }
223   
224    /**
225    * DOCUMENT ME!
226    *
227    * @return DOCUMENT ME!
228    */
 
229  147 toggle public String getAStr2()
230    {
231  147 return astr2;
232    }
233   
234    /**
235    * DOCUMENT ME!
236    *
237    * @return DOCUMENT ME!
238    */
 
239  2 toggle public int[] getASeq1()
240    {
241  2 return aseq1;
242    }
243   
244    /**
245    * DOCUMENT ME!
246    *
247    * @return DOCUMENT ME!
248    */
 
249  0 toggle public int[] getASeq2()
250    {
251  0 return aseq2;
252    }
253   
254    /**
255    *
256    * @return aligned instance of Seq 1
257    */
 
258  144 toggle public SequenceI getAlignedSeq1()
259    {
260  144 SequenceI alSeq1 = new Sequence(s1.getName(), getAStr1());
261  144 alSeq1.setStart(s1.getStart() + getSeq1Start() - 1);
262  144 alSeq1.setEnd(s1.getStart() + getSeq1End() - 1);
263  144 alSeq1.setDatasetSequence(
264  144 s1.getDatasetSequence() == null ? s1 : s1.getDatasetSequence());
265  144 return alSeq1;
266    }
267   
268    /**
269    *
270    * @return aligned instance of Seq 2
271    */
 
272  144 toggle public SequenceI getAlignedSeq2()
273    {
274  144 SequenceI alSeq2 = new Sequence(s2.getName(), getAStr2());
275  144 alSeq2.setStart(s2.getStart() + getSeq2Start() - 1);
276  144 alSeq2.setEnd(s2.getStart() + getSeq2End() - 1);
277  144 alSeq2.setDatasetSequence(
278  144 s2.getDatasetSequence() == null ? s2 : s2.getDatasetSequence());
279  144 return alSeq2;
280    }
281   
282    /**
283    * Construct score matrix for sequences with standard DNA or PEPTIDE matrix
284    *
285    * @param s1
286    * - sequence 1
287    * @param string1
288    * - string to use for s1
289    * @param s2
290    * - sequence 2
291    * @param string2
292    * - string to use for s2
293    * @param type
294    * DNA or PEPTIDE
295    */
 
296  181 toggle public void seqInit(SequenceI s1, String string1, SequenceI s2,
297    String string2, String type)
298    {
299  181 this.s1 = s1;
300  181 this.s2 = s2;
301  181 setDefaultParams(type);
302  181 seqInit(string1, string2);
303    }
304   
305    /**
306    * construct score matrix for string1 and string2 (after removing any existing
307    * gaps
308    *
309    * @param string1
310    * @param string2
311    */
 
312  181 toggle private void seqInit(String string1, String string2)
313    {
314  181 s1str = extractGaps(jalview.util.Comparison.GapChars, string1);
315  181 s2str = extractGaps(jalview.util.Comparison.GapChars, string2);
316   
317  181 if (s1str.length() == 0 || s2str.length() == 0)
318    {
319  0 output.append(
320  0 "ALL GAPS: " + (s1str.length() == 0 ? s1.getName() : " ")
321  0 + (s2str.length() == 0 ? s2.getName() : ""));
322  0 return;
323    }
324   
325  181 score = new float[s1str.length()][s2str.length()];
326   
327  181 E = new float[s1str.length()][s2str.length()];
328   
329  181 F = new float[s1str.length()][s2str.length()];
330  181 traceback = new int[s1str.length()][s2str.length()];
331   
332  181 seq1 = indexEncode(s1str);
333   
334  181 seq2 = indexEncode(s2str);
335    }
336   
 
337  181 toggle private void setDefaultParams(String moleculeType)
338    {
339  181 if (!PEP.equals(moleculeType) && !DNA.equals(moleculeType))
340    {
341  0 output.append("Wrong type = dna or pep only");
342  0 throw new Error(MessageManager
343    .formatMessage("error.unknown_type_dna_or_pep", new String[]
344    { moleculeType }));
345    }
346   
347  181 type = moleculeType;
348  181 scoreMatrix = ScoreModels.getInstance()
349    .getDefaultModel(PEP.equals(type));
350    }
351   
352    /**
353    * DOCUMENT ME!
354    */
 
355  179 toggle public void traceAlignment()
356    {
357    // Find the maximum score along the rhs or bottom row
358  179 float max = -Float.MAX_VALUE;
359   
360  21911 for (int i = 0; i < seq1.length; i++)
361    {
362  21732 if (score[i][seq2.length - 1] > max)
363    {
364  12067 max = score[i][seq2.length - 1];
365  12067 maxi = i;
366  12067 maxj = seq2.length - 1;
367    }
368    }
369   
370  19878 for (int j = 0; j < seq2.length; j++)
371    {
372  19699 if (score[seq1.length - 1][j] > max)
373    {
374  16 max = score[seq1.length - 1][j];
375  16 maxi = seq1.length - 1;
376  16 maxj = j;
377    }
378    }
379   
380  179 int i = maxi;
381  179 int j = maxj;
382  179 int trace;
383  179 maxscore = score[i][j] / 10f;
384   
385  179 seq1end = maxi + 1;
386  179 seq2end = maxj + 1;
387   
388  179 aseq1 = new int[seq1.length + seq2.length];
389  179 aseq2 = new int[seq1.length + seq2.length];
390   
391  179 StringBuilder sb1 = new StringBuilder(aseq1.length);
392  179 StringBuilder sb2 = new StringBuilder(aseq2.length);
393   
394  179 count = (seq1.length + seq2.length) - 1;
395   
396  16601 while (i > 0 && j > 0)
397    {
398  16422 aseq1[count] = seq1[i];
399  16422 sb1.append(s1str.charAt(i));
400  16422 aseq2[count] = seq2[j];
401  16422 sb2.append(s2str.charAt(j));
402   
403  16422 trace = findTrace(i, j);
404   
405  16422 if (trace == 0)
406    {
407  16422 i--;
408  16422 j--;
409    }
410  0 else if (trace == 1)
411    {
412  0 j--;
413  0 aseq1[count] = GAP_INDEX;
414  0 sb1.replace(sb1.length() - 1, sb1.length(), "-");
415    }
416  0 else if (trace == -1)
417    {
418  0 i--;
419  0 aseq2[count] = GAP_INDEX;
420  0 sb2.replace(sb2.length() - 1, sb2.length(), "-");
421    }
422   
423  16422 count--;
424    }
425   
426  179 seq1start = i + 1;
427  179 seq2start = j + 1;
428   
429  179 if (aseq1[count] != GAP_INDEX)
430    {
431  179 aseq1[count] = seq1[i];
432  179 sb1.append(s1str.charAt(i));
433    }
434   
435  179 if (aseq2[count] != GAP_INDEX)
436    {
437  179 aseq2[count] = seq2[j];
438  179 sb2.append(s2str.charAt(j));
439    }
440   
441    /*
442    * we built the character strings backwards, so now
443    * reverse them to convert to sequence strings
444    */
445  179 astr1 = sb1.reverse().toString();
446  179 astr2 = sb2.reverse().toString();
447    }
448   
449    /**
450    * DOCUMENT ME!
451    */
 
452  142 toggle public void printAlignment(PrintStream os)
453    {
454    // TODO: Use original sequence characters rather than re-translated
455    // characters in output
456    // Find the biggest id length for formatting purposes
457  142 String s1id = getAlignedSeq1().getDisplayId(true);
458  142 String s2id = getAlignedSeq2().getDisplayId(true);
459  142 int nameLength = Math.max(s1id.length(), s2id.length());
460  142 if (nameLength > MAX_NAME_LENGTH)
461    {
462  0 int truncateBy = nameLength - MAX_NAME_LENGTH;
463  0 nameLength = MAX_NAME_LENGTH;
464    // JAL-527 - truncate the sequence ids
465  0 if (s1id.length() > nameLength)
466    {
467  0 int slashPos = s1id.lastIndexOf('/');
468  0 s1id = s1id.substring(0, slashPos - truncateBy)
469    + s1id.substring(slashPos);
470    }
471  0 if (s2id.length() > nameLength)
472    {
473  0 int slashPos = s2id.lastIndexOf('/');
474  0 s2id = s2id.substring(0, slashPos - truncateBy)
475    + s2id.substring(slashPos);
476    }
477    }
478  142 int len = 72 - nameLength - 1;
479  142 int nochunks = ((aseq1.length - count) / len)
480  142 + ((aseq1.length - count) % len > 0 ? 1 : 0);
481  142 float pid = 0f;
482   
483  142 output.append("Score = ").append(score[maxi][maxj]).append(NEWLINE);
484  142 output.append("Length of alignment = ")
485    .append(String.valueOf(aseq1.length - count)).append(NEWLINE);
486  142 output.append("Sequence ");
487  142 Format nameFormat = new Format("%" + nameLength + "s");
488  142 output.append(nameFormat.form(s1id));
489  142 output.append(" (Sequence length = ")
490    .append(String.valueOf(s1str.length())).append(")")
491    .append(NEWLINE);
492  142 output.append("Sequence ");
493  142 output.append(nameFormat.form(s2id));
494  142 output.append(" (Sequence length = ")
495    .append(String.valueOf(s2str.length())).append(")")
496    .append(NEWLINE).append(NEWLINE);
497   
498  142 ScoreMatrix pam250 = ScoreModels.getInstance().getPam250();
499   
500  434 for (int j = 0; j < nochunks; j++)
501    {
502    // Print the first aligned sequence
503  292 output.append(nameFormat.form(s1id)).append(" ");
504   
505  16945 for (int i = 0; i < len; i++)
506    {
507  16653 if ((i + (j * len)) < astr1.length())
508    {
509  13547 output.append(astr1.charAt(i + (j * len)));
510    }
511    }
512   
513  292 output.append(NEWLINE);
514  292 output.append(nameFormat.form(" ")).append(" ");
515   
516    /*
517    * Print out the match symbols:
518    * | for exact match (ignoring case)
519    * . if PAM250 score is positive
520    * else a space
521    */
522  16945 for (int i = 0; i < len; i++)
523    {
524  16653 if ((i + (j * len)) < astr1.length())
525    {
526  13547 char c1 = astr1.charAt(i + (j * len));
527  13547 char c2 = astr2.charAt(i + (j * len));
528  13547 boolean sameChar = Comparison.isSameResidue(c1, c2, false);
529  13547 if (sameChar && !Comparison.isGap(c1))
530    {
531  13479 pid++;
532  13479 output.append("|");
533    }
534  68 else if (PEP.equals(type))
535    {
536  68 if (pam250.getPairwiseScore(c1, c2) > 0)
537    {
538  2 output.append(".");
539    }
540    else
541    {
542  66 output.append(" ");
543    }
544    }
545    else
546    {
547  0 output.append(" ");
548    }
549    }
550    }
551   
552    // Now print the second aligned sequence
553  292 output = output.append(NEWLINE);
554  292 output = output.append(nameFormat.form(s2id)).append(" ");
555   
556  16945 for (int i = 0; i < len; i++)
557    {
558  16653 if ((i + (j * len)) < astr2.length())
559    {
560  13547 output.append(astr2.charAt(i + (j * len)));
561    }
562    }
563   
564  292 output.append(NEWLINE).append(NEWLINE);
565    }
566   
567  142 pid = pid / (aseq1.length - count) * 100;
568  142 output.append(new Format("Percentage ID = %3.2f\n").form(pid));
569  142 output.append(NEWLINE);
570  142 try
571    {
572  142 os.print(output.toString());
573    } catch (Exception ex)
574    {
575    }
576    }
577   
578    /**
579    * DOCUMENT ME!
580    *
581    * @param i
582    * DOCUMENT ME!
583    * @param j
584    * DOCUMENT ME!
585    *
586    * @return DOCUMENT ME!
587    */
 
588  3093018 toggle public int findTrace(int i, int j)
589    {
590  3093018 int t = 0;
591  3093018 float pairwiseScore = scoreMatrix.getPairwiseScore(s1str.charAt(i),
592    s2str.charAt(j));
593  3093018 float max = score[i - 1][j - 1] + (pairwiseScore * 10);
594   
595  3093018 if (F[i][j] > max)
596    {
597  936537 max = F[i][j];
598  936537 t = -1;
599    }
600  2156481 else if (F[i][j] == max)
601    {
602  120320 if (prev == -1)
603    {
604  67095 max = F[i][j];
605  67095 t = -1;
606    }
607    }
608   
609  3093018 if (E[i][j] >= max)
610    {
611  1039755 max = E[i][j];
612  1039755 t = 1;
613    }
614  2053263 else if (E[i][j] == max)
615    {
616  0 if (prev == 1)
617    {
618  0 max = E[i][j];
619  0 t = 1;
620    }
621    }
622   
623  3093018 prev = t;
624   
625  3093018 return t;
626    }
627   
628    /**
629    * DOCUMENT ME!
630    */
 
631  179 toggle public void calcScoreMatrix()
632    {
633  179 int n = seq1.length;
634  179 int m = seq2.length;
635   
636    // top left hand element
637  179 score[0][0] = scoreMatrix.getPairwiseScore(s1str.charAt(0),
638    s2str.charAt(0)) * 10;
639  179 E[0][0] = -GAP_EXTEND_COST;
640  179 F[0][0] = 0;
641   
642    // Calculate the top row first
643  19699 for (int j = 1; j < m; j++)
644    {
645    // What should these values be? 0 maybe
646  19520 E[0][j] = max(score[0][j - 1] - GAP_OPEN_COST, E[0][j - 1] - GAP_EXTEND_COST);
647  19520 F[0][j] = -GAP_EXTEND_COST;
648   
649  19520 float pairwiseScore = scoreMatrix.getPairwiseScore(s1str.charAt(0),
650    s2str.charAt(j));
651  19520 score[0][j] = max(pairwiseScore * 10, -GAP_OPEN_COST, -GAP_EXTEND_COST);
652   
653  19520 traceback[0][j] = 1;
654    }
655   
656    // Now do the left hand column
657  21732 for (int i = 1; i < n; i++)
658    {
659  21553 E[i][0] = -GAP_OPEN_COST;
660  21553 F[i][0] = max(score[i - 1][0] - GAP_OPEN_COST, F[i - 1][0] - GAP_EXTEND_COST);
661   
662  21553 float pairwiseScore = scoreMatrix.getPairwiseScore(s1str.charAt(i),
663    s2str.charAt(0));
664  21553 score[i][0] = max(pairwiseScore * 10, E[i][0], F[i][0]);
665  21553 traceback[i][0] = -1;
666    }
667   
668    // Now do all the other rows
669  21732 for (int i = 1; i < n; i++)
670    {
671  3098149 for (int j = 1; j < m; j++)
672    {
673  3076596 E[i][j] = max(score[i][j - 1] - GAP_OPEN_COST, E[i][j - 1] - GAP_EXTEND_COST);
674  3076596 F[i][j] = max(score[i - 1][j] - GAP_OPEN_COST, F[i - 1][j] - GAP_EXTEND_COST);
675   
676  3076596 float pairwiseScore = scoreMatrix.getPairwiseScore(s1str.charAt(i),
677    s2str.charAt(j));
678  3076596 score[i][j] = max(score[i - 1][j - 1] + (pairwiseScore * 10),
679    E[i][j], F[i][j]);
680  3076596 traceback[i][j] = findTrace(i, j);
681    }
682    }
683    }
684   
685    /**
686    * Returns the given sequence with all of the given gap characters removed.
687    *
688    * @param gapChars
689    * a string of characters to be treated as gaps
690    * @param seq
691    * the input sequence
692    *
693    * @return
694    */
 
695  4198 toggle public static String extractGaps(String gapChars, String seq)
696    {
697  4198 if (gapChars == null || seq == null)
698    {
699  6 return null;
700    }
701  4192 StringTokenizer str = new StringTokenizer(seq, gapChars);
702  4192 StringBuilder newString = new StringBuilder(seq.length());
703   
704  24769 while (str.hasMoreTokens())
705    {
706  20577 newString.append(str.nextToken());
707    }
708   
709  4192 return newString.toString();
710    }
711   
712    /**
713    * DOCUMENT ME!
714    *
715    * @param f1
716    * DOCUMENT ME!
717    * @param f2
718    * DOCUMENT ME!
719    * @param f3
720    * DOCUMENT ME!
721    *
722    * @return DOCUMENT ME!
723    */
 
724  3117669 toggle private static float max(float f1, float f2, float f3)
725    {
726  3117669 float max = f1;
727   
728  3117669 if (f2 > f1)
729    {
730  941269 max = f2;
731    }
732   
733  3117669 if (f3 > max)
734    {
735  928010 max = f3;
736    }
737   
738  3117669 return max;
739    }
740   
741    /**
742    * DOCUMENT ME!
743    *
744    * @param f1
745    * DOCUMENT ME!
746    * @param f2
747    * DOCUMENT ME!
748    *
749    * @return DOCUMENT ME!
750    */
 
751  6194265 toggle private static float max(float f1, float f2)
752    {
753  6194265 float max = f1;
754   
755  6194265 if (f2 > f1)
756    {
757  3314432 max = f2;
758    }
759   
760  6194265 return max;
761    }
762   
763    /**
764    * Converts the character string to an array of integers which are the
765    * corresponding indices to the characters in the score matrix
766    *
767    * @param s
768    *
769    * @return
770    */
 
771  364 toggle int[] indexEncode(String s)
772    {
773  364 int[] encoded = new int[s.length()];
774   
775  41848 for (int i = 0; i < s.length(); i++)
776    {
777  41484 char c = s.charAt(i);
778  41484 encoded[i] = scoreMatrix.getMatrixIndex(c);
779    }
780   
781  364 return encoded;
782    }
783   
784    /**
785    * DOCUMENT ME!
786    *
787    * @param g
788    * DOCUMENT ME!
789    * @param mat
790    * DOCUMENT ME!
791    * @param n
792    * DOCUMENT ME!
793    * @param m
794    * DOCUMENT ME!
795    * @param psize
796    * DOCUMENT ME!
797    */
 
798  0 toggle public static void displayMatrix(Graphics g, int[][] mat, int n, int m,
799    int psize)
800    {
801    // TODO method doesn't seem to be referenced anywhere delete??
802  0 int max = -1000;
803  0 int min = 1000;
804   
805  0 for (int i = 0; i < n; i++)
806    {
807  0 for (int j = 0; j < m; j++)
808    {
809  0 if (mat[i][j] >= max)
810    {
811  0 max = mat[i][j];
812    }
813   
814  0 if (mat[i][j] <= min)
815    {
816  0 min = mat[i][j];
817    }
818    }
819    }
820   
821  0 System.out.println(max + " " + min);
822   
823  0 for (int i = 0; i < n; i++)
824    {
825  0 for (int j = 0; j < m; j++)
826    {
827  0 int x = psize * i;
828  0 int y = psize * j;
829   
830    // System.out.println(mat[i][j]);
831  0 float score = (float) (mat[i][j] - min) / (float) (max - min);
832  0 g.setColor(new Color(score, 0, 0));
833  0 g.fillRect(x, y, psize, psize);
834   
835    // System.out.println(x + " " + y + " " + score);
836    }
837    }
838    }
839   
840    /**
841    * Compute a globally optimal needleman and wunsch alignment between two
842    * sequences
843    *
844    * @param s1
845    * @param s2
846    * @param type
847    * AlignSeq.DNA or AlignSeq.PEP
848    */
 
849  177 toggle public static AlignSeq doGlobalNWAlignment(SequenceI s1, SequenceI s2,
850    String type)
851    {
852  177 AlignSeq as = new AlignSeq(s1, s2, type);
853   
854  177 as.calcScoreMatrix();
855  177 as.traceAlignment();
856  177 return as;
857    }
858   
859    /**
860    *
861    * @return mapping from positions in S1 to corresponding positions in S2
862    */
 
863  169 toggle public jalview.datamodel.Mapping getMappingFromS1(boolean allowmismatch)
864    {
865  169 ArrayList<Integer> as1 = new ArrayList<Integer>(),
866    as2 = new ArrayList<Integer>();
867  169 int pdbpos = s2.getStart() + getSeq2Start() - 2;
868  169 int alignpos = s1.getStart() + getSeq1Start() - 2;
869  169 int lp2 = pdbpos - 3, lp1 = alignpos - 3;
870  169 boolean lastmatch = false;
871    // and now trace the alignment onto the atom set.
872  16912 for (int i = 0; i < astr1.length(); i++)
873    {
874  16743 char c1 = astr1.charAt(i), c2 = astr2.charAt(i);
875  16743 if (c1 != '-')
876    {
877  16743 alignpos++;
878    }
879   
880  16743 if (c2 != '-')
881    {
882  16743 pdbpos++;
883    }
884   
885  16743 if (allowmismatch || c1 == c2)
886    {
887    // extend mapping interval
888  16681 if (lp1 + 1 != alignpos || lp2 + 1 != pdbpos)
889    {
890  230 as1.add(Integer.valueOf(alignpos));
891  230 as2.add(Integer.valueOf(pdbpos));
892    }
893  16681 lastmatch = true;
894  16681 lp1 = alignpos;
895  16681 lp2 = pdbpos;
896    }
897    else
898    {
899    // extend mapping interval
900  62 if (lastmatch)
901    {
902  61 as1.add(Integer.valueOf(lp1));
903  61 as2.add(Integer.valueOf(lp2));
904    }
905  62 lastmatch = false;
906    }
907    }
908    // construct range pairs
909   
910  169 int[] mapseq1 = new int[as1.size() + (lastmatch ? 1 : 0)],
911  169 mapseq2 = new int[as2.size() + (lastmatch ? 1 : 0)];
912  169 int i = 0;
913  169 for (Integer ip : as1)
914    {
915  291 mapseq1[i++] = ip;
916    }
917  169 ;
918  169 i = 0;
919  169 for (Integer ip : as2)
920    {
921  291 mapseq2[i++] = ip;
922    }
923  169 ;
924  169 if (lastmatch)
925    {
926  169 mapseq1[mapseq1.length - 1] = alignpos;
927  169 mapseq2[mapseq2.length - 1] = pdbpos;
928    }
929  169 MapList map = new MapList(mapseq1, mapseq2, 1, 1);
930   
931  169 jalview.datamodel.Mapping mapping = new Mapping(map);
932  169 mapping.setTo(s2);
933  169 return mapping;
934    }
935   
936    /**
937    * matches ochains against al and populates seqs with the best match between
938    * each ochain and the set in al
939    *
940    * @param ochains
941    * @param al
942    * @param dnaOrProtein
943    * @param removeOldAnnots
944    * when true, old annotation is cleared before new annotation
945    * transferred
946    * @return List<List<SequenceI> originals, List<SequenceI> replacement,
947    * List<AlignSeq> alignment between each>
948    */
 
949  2 toggle public static List<List<? extends Object>> replaceMatchingSeqsWith(
950    List<SequenceI> seqs, List<AlignmentAnnotation> annotations,
951    List<SequenceI> ochains, AlignmentI al, String dnaOrProtein,
952    boolean removeOldAnnots)
953    {
954  2 List<SequenceI> orig = new ArrayList<SequenceI>(),
955    repl = new ArrayList<SequenceI>();
956  2 List<AlignSeq> aligs = new ArrayList<AlignSeq>();
957  2 if (al != null && al.getHeight() > 0)
958    {
959  2 ArrayList<SequenceI> matches = new ArrayList<SequenceI>();
960  2 ArrayList<AlignSeq> aligns = new ArrayList<AlignSeq>();
961   
962  2 for (SequenceI sq : ochains)
963    {
964  8 SequenceI bestm = null;
965  8 AlignSeq bestaseq = null;
966  8 float bestscore = 0;
967  8 for (SequenceI msq : al.getSequences())
968    {
969  20 AlignSeq aseq = doGlobalNWAlignment(msq, sq, dnaOrProtein);
970  20 if (bestm == null || aseq.getMaxScore() > bestscore)
971    {
972  8 bestscore = aseq.getMaxScore();
973  8 bestaseq = aseq;
974  8 bestm = msq;
975    }
976    }
977    // System.out.println("Best Score for " + (matches.size() + 1) + " :"
978    // + bestscore);
979  8 matches.add(bestm);
980  8 aligns.add(bestaseq);
981  8 al.deleteSequence(bestm);
982    }
983  10 for (int p = 0, pSize = seqs.size(); p < pSize; p++)
984    {
985  8 SequenceI sq, sp = seqs.get(p);
986  8 int q;
987  ? if ((q = ochains.indexOf(sp)) > -1)
988    {
989  8 seqs.set(p, sq = matches.get(q));
990  8 orig.add(sp);
991  8 repl.add(sq);
992  8 sq.setName(sp.getName());
993  8 sq.setDescription(sp.getDescription());
994  8 Mapping sp2sq;
995  8 sq.transferAnnotation(sp,
996    sp2sq = aligns.get(q).getMappingFromS1(false));
997  8 aligs.add(aligns.get(q));
998  8 int inspos = -1;
999  36 for (int ap = 0; ap < annotations.size();)
1000    {
1001  28 if (annotations.get(ap).sequenceRef == sp)
1002    {
1003  4 if (inspos == -1)
1004    {
1005  4 inspos = ap;
1006    }
1007  4 if (removeOldAnnots)
1008    {
1009  0 annotations.remove(ap);
1010    }
1011    else
1012    {
1013  4 AlignmentAnnotation alan = annotations.remove(ap);
1014  4 alan.liftOver(sq, sp2sq);
1015  4 alan.setSequenceRef(sq);
1016  4 sq.addAlignmentAnnotation(alan);
1017    }
1018    }
1019    else
1020    {
1021  24 ap++;
1022    }
1023    }
1024  8 if (sq.getAnnotation() != null && sq.getAnnotation().length > 0)
1025    {
1026  8 annotations.addAll(inspos == -1 ? annotations.size() : inspos,
1027    Arrays.asList(sq.getAnnotation()));
1028    }
1029    }
1030    }
1031    }
1032  2 return Arrays.asList(orig, repl, aligs);
1033    }
1034   
1035    /**
1036    * compute the PID vector used by the redundancy filter.
1037    *
1038    * @param originalSequences
1039    * - sequences in alignment that are to filtered
1040    * @param omitHidden
1041    * - null or strings to be analysed (typically, visible portion of
1042    * each sequence in alignment)
1043    * @param start
1044    * - first column in window for calculation
1045    * @param end
1046    * - last column in window for calculation
1047    * @param ungapped
1048    * - if true then use ungapped sequence to compute PID
1049    * @return vector containing maximum PID for i-th sequence and any sequences
1050    * longer than that seuqence
1051    */
 
1052  0 toggle public static float[] computeRedundancyMatrix(
1053    SequenceI[] originalSequences, String[] omitHidden, int start,
1054    int end, boolean ungapped)
1055    {
1056  0 int height = originalSequences.length;
1057  0 float[] redundancy = new float[height];
1058  0 int[] lngth = new int[height];
1059  0 for (int i = 0; i < height; i++)
1060    {
1061  0 redundancy[i] = 0f;
1062  0 lngth[i] = -1;
1063    }
1064   
1065    // long start = System.currentTimeMillis();
1066   
1067  0 SimilarityParams pidParams = new SimilarityParams(true, true, true,
1068    true);
1069  0 float pid;
1070  0 String seqi, seqj;
1071  0 for (int i = 0; i < height; i++)
1072    {
1073   
1074  0 for (int j = 0; j < i; j++)
1075    {
1076  0 if (i == j)
1077    {
1078  0 continue;
1079    }
1080   
1081  0 if (omitHidden == null)
1082    {
1083  0 seqi = originalSequences[i].getSequenceAsString(start, end);
1084  0 seqj = originalSequences[j].getSequenceAsString(start, end);
1085    }
1086    else
1087    {
1088  0 seqi = omitHidden[i];
1089  0 seqj = omitHidden[j];
1090    }
1091  0 if (lngth[i] == -1)
1092    {
1093  0 String ug = AlignSeq.extractGaps(Comparison.GapChars, seqi);
1094  0 lngth[i] = ug.length();
1095  0 if (ungapped)
1096    {
1097  0 seqi = ug;
1098    }
1099    }
1100  0 if (lngth[j] == -1)
1101    {
1102  0 String ug = AlignSeq.extractGaps(Comparison.GapChars, seqj);
1103  0 lngth[j] = ug.length();
1104  0 if (ungapped)
1105    {
1106  0 seqj = ug;
1107    }
1108    }
1109  0 pid = (float) PIDModel.computePID(seqi, seqj, pidParams);
1110   
1111    // use real sequence length rather than string length
1112  0 if (lngth[j] < lngth[i])
1113    {
1114  0 redundancy[j] = Math.max(pid, redundancy[j]);
1115    }
1116    else
1117    {
1118  0 redundancy[i] = Math.max(pid, redundancy[i]);
1119    }
1120   
1121    }
1122    }
1123  0 return redundancy;
1124    }
1125    }