src/jalview/datamodel/SeqCigar.java

   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.datamodel;
  22
  23 import jalview.analysis.AlignSeq;
  24 import jalview.analysis.SeqsetUtils;
  25 import jalview.util.MessageManager;
  26 import jalview.util.ShiftList;
  27
  28 import java.util.Enumeration;
  29 import java.util.Hashtable;
  30
  31 public class SeqCigar extends CigarSimple
  32 {
  33   /**
  34    * start(inclusive) and end(exclusive) of subsequence on refseq
  35    */
  36   private int start, end;
  37
  38   private SequenceI refseq = null;
  39
  40   private Hashtable seqProps;
  41
  42   /**
  43    * Reference dataset sequence for the cigar string
  44    *
  45    * @return SequenceI
  46    */
  47   public SequenceI getRefSeq()
  48   {
  49     return refseq;
  50   }
  51
  52   /**
  53    *
  54    * @return int start index of cigar ops on refSeq
  55    */
  56   public int getStart()
  57   {
  58     return start;
  59   }
  60
  61   /**
  62    *
  63    * @return int end index (exclusive) of cigar ops on refSeq
  64    */
  65   public int getEnd()
  66   {
  67     return end;
  68   }
  69
  70   /**
  71   /**
  72    * Returns sequence as a string with cigar operations applied to it
  73    *
  74    * @return String
  75    */
  76   @Override
  77   public String getSequenceString(char GapChar)
  78   {
  79     return (length == 0) ? "" : (String) getSequenceAndDeletions(
  80             refseq.getSequenceAsString(start, end), GapChar)[0];
  81   }
  82
  83   /**
  84    * recreates a gapped and edited version of RefSeq or null for an empty cigar
  85    * string
  86    *
  87    * @return SequenceI
  88    */
  89   public SequenceI getSeq(char GapChar)
  90   {
  91     Sequence seq;
  92     if (refseq == null || length == 0)
  93     {
  94       return null;
  95     }
  96     Object[] edit_result = getSequenceAndDeletions(
  97             refseq.getSequenceAsString(start, end), GapChar);
  98     if (edit_result == null)
  99     {
 100       throw new Error(
 101               MessageManager
 102                       .getString("error.implementation_error_unexpected_null_from_get_sequence_and_deletions"));
 103     }
 104     int bounds[] = (int[]) edit_result[1];
 105     seq = new Sequence(refseq.getName(), (String) edit_result[0],
 106             refseq.getStart() + start + bounds[0], refseq.getStart()
 107                     + start + ((bounds[2] == 0) ? -1 : bounds[2]));
 108     seq.setDescription(refseq.getDescription());
 109     int sstart = seq.getStart(), send = seq.getEnd();
 110     // seq.checkValidRange(); probably not needed
 111     // recover local properties if present
 112     if (seqProps != null)
 113     {
 114       // this recovers dataset sequence reference as well as local features,
 115       // names, start/end settings.
 116       SeqsetUtils.SeqCharacterUnhash(seq, seqProps);
 117     }
 118     // ensure dataset sequence is up to date from local reference
 119     seq.setDatasetSequence(refseq);
 120     seq.setStart(sstart);
 121     seq.setEnd(send);
 122     return seq;
 123   }
 124
 125   /*
 126    * We don't allow this - refseq is given at construction time only public void
 127    * setSeq(SequenceI seq) { this.seq = seq; }
 128    */
 129   /**
 130    * internal constructor - sets seq to a gapless sequence derived from seq and
 131    * prepends any 'D' operations needed to get to the first residue of seq.
 132    *
 133    * @param seq
 134    *          SequenceI
 135    * @param initialDeletion
 136    *          true to mark initial dataset sequence residues as deleted in
 137    *          subsequence
 138    * @param _s
 139    *          index of first position in seq
 140    * @param _e
 141    *          index after last position in (possibly gapped) seq
 142    * @return true if gaps are present in seq
 143    */
 144   private boolean _setSeq(SequenceI seq, boolean initialDeletion, int _s,
 145           int _e)
 146   {
 147     boolean hasgaps = false;
 148     if (seq == null)
 149     {
 150       throw new Error(
 151               MessageManager
 152                       .getString("error.implementation_error_set_seq_null"));
 153     }
 154     if (_s < 0)
 155     {
 156       throw new Error(MessageManager.formatMessage(
 157               "error.implementation_error_s", new String[] { Integer
 158                       .valueOf(_s).toString() }));
 159     }
 160     String seq_string = seq.getSequenceAsString();
 161     if (_e == 0 || _e < _s || _e > seq_string.length())
 162     {
 163       _e = seq_string.length();
 164     }
 165     // resolve start and end positions relative to ungapped reference sequence
 166     start = seq.findPosition(_s) - seq.getStart();
 167     end = seq.findPosition(_e) - seq.getStart();
 168     int l_ungapped = end - start;
 169     // Find correct sequence to reference and correct start and end - if
 170     // necessary
 171     SequenceI ds = seq.getDatasetSequence();
 172     if (ds == null)
 173     {
 174       // make a new dataset sequence
 175       String ungapped = AlignSeq.extractGaps(
 176               jalview.util.Comparison.GapChars, new String(seq_string));
 177       l_ungapped = ungapped.length();
 178       // check that we haven't just duplicated an ungapped sequence.
 179       if (l_ungapped == seq.getLength())
 180       {
 181         ds = seq;
 182       }
 183       else
 184       {
 185         ds = new Sequence(seq.getName(), ungapped, seq.getStart(),
 186                 seq.getStart() + ungapped.length() - 1);
 187         // JBPNote: this would be consistent but may not be useful
 188         // seq.setDatasetSequence(ds);
 189       }
 190     }
 191     // add in offset between seq and the dataset sequence
 192     if (ds.getStart() < seq.getStart())
 193     {
 194       int offset = seq.getStart() - ds.getStart();
 195       if (initialDeletion)
 196       {
 197         // absolute cigar string
 198         addDeleted(_s + offset);
 199         start = 0;
 200         end += offset;
 201       }
 202       else
 203       {
 204         // normal behaviour - just mark start and end subsequence
 205         start += offset;
 206         end += offset;
 207
 208       }
 209
 210     }
 211
 212     // any gaps to process ?
 213     if (l_ungapped != (_e - _s))
 214     {
 215       hasgaps = true;
 216     }
 217
 218     refseq = ds;
 219     // copy over local properties for the sequence instance of the refseq
 220     seqProps = SeqsetUtils.SeqCharacterHash(seq);
 221     // Check offsets
 222     if (end > ds.getLength())
 223     {
 224       throw new Error(
 225               MessageManager
 226                       .getString("error.implementation_error_seqcigar_possible"));
 227       // end = ds.getLength();
 228     }
 229
 230     return hasgaps;
 231   }
 232
 233   /**
 234    * directly initialise a cigar object with a sequence of range, operation
 235    * pairs and a sequence to apply it to. operation and range should be relative
 236    * to the seq.getStart()'th residue of the dataset seq resolved from seq.
 237    *
 238    * @param seq
 239    *          SequenceI
 240    * @param operation
 241    *          char[]
 242    * @param range
 243    *          int[]
 244    */
 245   public SeqCigar(SequenceI seq, char operation[], int range[])
 246   {
 247     super();
 248     if (seq == null)
 249     {
 250       throw new Error(
 251               MessageManager.getString("error.implmentation_bug_seq_null"));
 252     }
 253     if (operation.length != range.length)
 254     {
 255       throw new Error(
 256               MessageManager
 257                       .getString("error.implementation_bug_cigar_operation_list_range_list"));
 258     }
 259
 260     if (operation != null)
 261     {
 262       this.operation = new char[operation.length + _inc_length];
 263       this.range = new int[operation.length + _inc_length];
 264
 265       if (_setSeq(seq, false, 0, 0))
 266       {
 267         throw new Error(
 268                 MessageManager
 269                         .getString("error.not_yet_implemented_cigar_object_from_cigar_string"));
 270       }
 271       for (int i = this.length, j = 0; j < operation.length; i++, j++)
 272       {
 273         char op = operation[j];
 274         if (op != M && op != I && op != D)
 275         {
 276           throw new Error(MessageManager.formatMessage(
 277                   "error.implementation_bug_cigar_operation", new String[] {
 278                       Integer.valueOf(j).toString(),
 279                       Integer.valueOf(op).toString(),
 280                       Integer.valueOf(M).toString(),
 281                       Integer.valueOf(I).toString(),
 282                       Integer.valueOf(D).toString() }));
 283         }
 284         this.operation[i] = op;
 285         this.range[i] = range[j];
 286       }
 287       this.length += operation.length;
 288     }
 289     else
 290     {
 291       this.operation = null;
 292       this.range = null;
 293       this.length = 0;
 294       if (_setSeq(seq, false, 0, 0))
 295       {
 296         throw new Error(
 297                 MessageManager
 298                         .getString("error.not_yet_implemented_cigar_object_from_cigar_string"));
 299       }
 300     }
 301   }
 302
 303   /**
 304    * add range matched residues to cigar string
 305    *
 306    * @param range
 307    *          int
 308    */
 309   public void addMatch(int range)
 310   {
 311     this.addOperation(M, range);
 312   }
 313
 314   /**
 315    * Adds insertion and match operations based on seq to the cigar up to the
 316    * endpos column of seq.
 317    *
 318    * @param cigar
 319    *          CigarBase
 320    * @param seq
 321    *          SequenceI
 322    * @param startpos
 323    *          int
 324    * @param endpos
 325    *          int
 326    * @param initialDeletions
 327    *          if true then initial deletions will be added from start of seq to
 328    *          startpos
 329    */
 330   protected static void addSequenceOps(CigarBase cigar, SequenceI seq,
 331           int startpos, int endpos, boolean initialDeletions)
 332   {
 333     char op = '\0';
 334     int range = 0;
 335     int p = 0, res = seq.getLength();
 336
 337     if (!initialDeletions)
 338     {
 339       p = startpos;
 340     }
 341
 342     while (p <= endpos)
 343     {
 344       boolean isGap = (p < res) ? jalview.util.Comparison.isGap(seq
 345               .getCharAt(p)) : true;
 346       if ((startpos <= p) && (p <= endpos))
 347       {
 348         if (isGap)
 349         {
 350           if (range > 0 && op != I)
 351           {
 352             cigar.addOperation(op, range);
 353             range = 0;
 354           }
 355           op = I;
 356           range++;
 357         }
 358         else
 359         {
 360           if (range > 0 && op != M)
 361           {
 362             cigar.addOperation(op, range);
 363             range = 0;
 364           }
 365           op = M;
 366           range++;
 367         }
 368       }
 369       else
 370       {
 371         if (!isGap)
 372         {
 373           if (range > 0 && op != D)
 374           {
 375             cigar.addOperation(op, range);
 376             range = 0;
 377           }
 378           op = D;
 379           range++;
 380         }
 381         else
 382         {
 383           // do nothing - insertions are not made in flanking regions
 384         }
 385       }
 386       p++;
 387     }
 388     if (range > 0)
 389     {
 390       cigar.addOperation(op, range);
 391     }
 392   }
 393
 394   /**
 395    * create a cigar string for given sequence
 396    *
 397    * @param seq
 398    *          SequenceI
 399    */
 400   public SeqCigar(SequenceI seq)
 401   {
 402     super();
 403     if (seq == null)
 404     {
 405       throw new Error(
 406               MessageManager
 407                       .getString("error.implementation_error_for_new_cigar"));
 408     }
 409     _setSeq(seq, false, 0, 0);
 410     // there is still work to do
 411     addSequenceOps(this, seq, 0, seq.getLength() - 1, false);
 412   }
 413
 414   /**
 415    * Create Cigar from a range of gaps and residues on a sequence object
 416    *
 417    * @param seq
 418    *          SequenceI
 419    * @param start
 420    *          int - first column in range
 421    * @param end
 422    *          int - last column in range
 423    */
 424   public SeqCigar(SequenceI seq, int start, int end)
 425   {
 426     super();
 427     if (seq == null)
 428     {
 429       throw new Error(
 430               MessageManager
 431                       .getString("error.implementation_error_for_new_cigar"));
 432     }
 433     _setSeq(seq, false, start, end + 1);
 434     // there is still work to do
 435     addSequenceOps(this, seq, start, end, false);
 436   }
 437
 438   /**
 439    * Create a cigar object from a cigar string like '[<I|D|M><range>]+' Will
 440    * fail if the given seq already contains gaps (JBPNote: future implementation
 441    * will fix)
 442    *
 443    * @param seq
 444    *          SequenceI object resolvable to a dataset sequence
 445    * @param cigarString
 446    *          String
 447    * @return Cigar
 448    */
 449   public static SeqCigar parseCigar(SequenceI seq, String cigarString)
 450           throws Exception
 451   {
 452     Object[] opsandrange = parseCigarString(cigarString);
 453     return new SeqCigar(seq, (char[]) opsandrange[0],
 454             (int[]) opsandrange[1]);
 455   }
 456
 457   /**
 458    * create an alignment from the given array of cigar sequences and gap
 459    * character, and marking the given segments as visible in the given
 460    * columselection.
 461    *
 462    * @param alseqs
 463    * @param gapCharacter
 464    * @param colsel
 465    *          - columnSelection where hidden regions are marked
 466    * @param segments
 467    *          - visible regions of alignment
 468    * @return SequenceI[]
 469    */
 470   public static SequenceI[] createAlignmentSequences(SeqCigar[] alseqs,
 471           char gapCharacter, ColumnSelection colsel, int[] segments)
 472   {
 473     SequenceI[] seqs = new SequenceI[alseqs.length];
 474     StringBuffer[] g_seqs = new StringBuffer[alseqs.length];
 475     String[] alseqs_string = new String[alseqs.length];
 476     Object[] gs_regions = new Object[alseqs.length];
 477     for (int i = 0; i < alseqs.length; i++)
 478     {
 479       alseqs_string[i] = alseqs[i].getRefSeq().getSequenceAsString(
 480               alseqs[i].start, alseqs[i].end);
 481       gs_regions[i] = alseqs[i].getSequenceAndDeletions(alseqs_string[i],
 482               gapCharacter); // gapped sequence, {start, start col, end.
 483       // endcol}, hidden regions {{start, end, col}})
 484       if (gs_regions[i] == null)
 485       {
 486         throw new Error(MessageManager.formatMessage(
 487                 "error.implementation_error_cigar_seq_no_operations",
 488                 new String[] { Integer.valueOf(i).toString() }));
 489       }
 490       g_seqs[i] = new StringBuffer((String) ((Object[]) gs_regions[i])[0]); // the
 491       // visible
 492       // gapped
 493       // sequence
 494     }
 495     // Now account for insertions. (well - deletions)
 496     // this is complicated because we must keep track of shifted positions in
 497     // each sequence
 498     ShiftList shifts = new ShiftList();
 499     for (int i = 0; i < alseqs.length; i++)
 500     {
 501       Object[] gs_region = ((Object[]) ((Object[]) gs_regions[i])[2]);
 502       if (gs_region != null)
 503
 504       {
 505         for (int hr = 0; hr < gs_region.length; hr++)
 506         {
 507           int[] region = (int[]) gs_region[hr];
 508           char[] insert = new char[region[1] - region[0] + 1];
 509           for (int s = 0; s < insert.length; s++)
 510           {
 511             insert[s] = gapCharacter;
 512           }
 513           int inspos = shifts.shift(region[2]); // resolve insertion position in
 514           // current alignment frame of
 515           // reference
 516           for (int s = 0; s < alseqs.length; s++)
 517           {
 518             if (s != i)
 519             {
 520               if (g_seqs[s].length() <= inspos)
 521               {
 522                 // prefix insertion with more gaps.
 523                 for (int l = inspos - g_seqs[s].length(); l > 0; l--)
 524                 {
 525                   g_seqs[s].append(gapCharacter); // to debug - use a diffferent
 526                   // gap character here
 527                 }
 528               }
 529               g_seqs[s].insert(inspos, insert);
 530             }
 531             else
 532             {
 533               g_seqs[s].insert(inspos,
 534                       alseqs_string[i].substring(region[0], region[1] + 1));
 535             }
 536           }
 537           shifts.addShift(region[2], insert.length); // update shift in
 538           // alignment frame of
 539           // reference
 540           if (segments == null)
 541           {
 542             // add a hidden column for this deletion
 543             colsel.hideColumns(inspos, inspos + insert.length - 1);
 544           }
 545         }
 546       }
 547     }
 548     for (int i = 0; i < alseqs.length; i++)
 549     {
 550       int[] bounds = ((int[]) ((Object[]) gs_regions[i])[1]);
 551       SequenceI ref = alseqs[i].getRefSeq();
 552       seqs[i] = new Sequence(ref.getName(), g_seqs[i].toString(),
 553               ref.getStart() + alseqs[i].start + bounds[0], ref.getStart()
 554                       + alseqs[i].start + (bounds[2] == 0 ? -1 : bounds[2]));
 555       seqs[i].setDatasetSequence(ref);
 556       seqs[i].setDescription(ref.getDescription());
 557     }
 558     if (segments != null)
 559     {
 560       for (int i = 0; i < segments.length; i += 3)
 561       {
 562         // int start=shifts.shift(segments[i]-1)+1;
 563         // int end=shifts.shift(segments[i]+segments[i+1]-1)-1;
 564         colsel.hideColumns(segments[i + 1], segments[i + 1]
 565                 + segments[i + 2] - 1);
 566       }
 567     }
 568     return seqs;
 569   }
 570
 571   /**
 572    * references to entities that this sequence cigar is associated with.
 573    */
 574   private Hashtable selGroups = null;
 575
 576   public void setGroupMembership(Object group)
 577   {
 578     if (selGroups == null)
 579     {
 580       selGroups = new Hashtable();
 581     }
 582     selGroups.put(group, new int[0]);
 583   }
 584
 585   /**
 586    * Test for and if present remove association to group.
 587    *
 588    * @param group
 589    * @return true if group was associated and it was removed
 590    */
 591   public boolean removeGroupMembership(Object group)
 592   {
 593     if (selGroups != null && selGroups.containsKey(group))
 594     {
 595       selGroups.remove(group);
 596       return true;
 597     }
 598     return false;
 599   }
 600
 601   /**
 602    * forget all associations for this sequence.
 603    */
 604   public void clearMemberships()
 605   {
 606     if (selGroups != null)
 607     {
 608       selGroups.clear();
 609     }
 610     selGroups = null;
 611   }
 612
 613   /**
 614    *
 615    * @return null or array of all associated entities
 616    */
 617   public Object[] getAllMemberships()
 618   {
 619     if (selGroups == null)
 620     {
 621       return null;
 622     }
 623     Object[] mmbs = new Object[selGroups.size()];
 624     Enumeration en = selGroups.keys();
 625     for (int i = 0; en.hasMoreElements(); i++)
 626     {
 627       mmbs[i] = en.nextElement();
 628     }
 629     return mmbs;
 630   }
 631
 632   /**
 633    * Test for group membership
 634    *
 635    * @param sgr
 636    *          - a selection group or some other object that may be associated
 637    *          with seqCigar
 638    * @return true if sgr is associated with this seqCigar
 639    */
 640   public boolean isMemberOf(Object sgr)
 641   {
 642     return (selGroups != null) && selGroups.get(sgr) != null;
 643   }
 644 }