开通博客赚积分
发布资源赚积分
分类

源码开发语言/平台
当前位置：首页 > 源码/资料 > 数值算法/人工智能 > 生物技术 > 查看源码
cav_alndisplay.cpp
资源名称：ncbi_cxx_win--Jul_08_2004.zip　[点击查看]
上传用户：yhdzpy8989
上传日期：2007-06-13
资源大小：13604k
文件大小：51k
源码类别：

生物技术

开发平台：

C/C++

cav_alndisplay.cpp：源码内容
  1. /*
  2.  * ===========================================================================
  3.  * PRODUCTION $Log: cav_alndisplay.cpp,v $
  4.  * PRODUCTION Revision 1000.3  2004/06/01 19:41:15  gouriano
  5.  * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.5
  6.  * PRODUCTION
  7.  * ===========================================================================
  8.  */
  10. /*  $Id: cav_alndisplay.cpp,v 1000.3 2004/06/01 19:41:15 gouriano Exp $
  11. * ===========================================================================
  12. *
  13. *                            PUBLIC DOMAIN NOTICE
  14. *               National Center for Biotechnology Information
  15. *
  16. *  This software/database is a "United States Government Work" under the
  17. *  terms of the United States Copyright Act.  It was written as part of
  18. *  the author's official duties as a United States Government employee and
  19. *  thus cannot be copyrighted.  This software/database is freely available
  20. *  to the public for use. The National Library of Medicine and the U.S.
  21. *  Government have not placed any restriction on its use or reproduction.
  22. *
  23. *  Although all reasonable efforts have been taken to ensure the accuracy
  24. *  and reliability of the software and data, the NLM and the U.S.
  25. *  Government do not and cannot warrant the performance or results that
  26. *  may be obtained by using this software or data. The NLM and the U.S.
  27. *  Government disclaim all warranties, express or implied, including
  28. *  warranties of performance, merchantability or fitness for any particular
  29. *  purpose.
  30. *
  31. *  Please cite the author in any work or product based on this material.
  32. *
  33. * ===========================================================================
  34. *
  35. * Authors:  Paul Thiessen
  36. *
  37. * File Description:
  38. *      Classes to hold alignment display
  39. *
  40. * ===========================================================================
  41. */
  43. #include <ncbi_pch.hpp>
  44. #include <corelib/ncbistl.hpp>
  45. #include <corelib/ncbiobj.hpp>
  47. #include <map>
  48. #include <deque>
  49. #include <memory>
  50. #include <iomanip>
  51. #include <math.h>
  53. #include <objtools/cddalignview/cav_alndisplay.hpp>
  54. #include <objtools/cddalignview/cav_seqset.hpp>
  55. #include <objtools/cddalignview/cav_alignset.hpp>
  56. #include <objtools/cddalignview/cddalignview.h>
  59. BEGIN_NCBI_SCOPE
  61. const double AlignmentDisplay::SHOW_IDENTITY = 100000.0;
  63. // HTML colors
  64. static const int nBlockColors = 2;
  65. static const string
  66.     bgColor("#FFFFFF"), blockBGColors[nBlockColors] = { "#FFFFE0", "#FFE8FF" },
  67.     rulerColor("#700777"), numColor("#229922"), featColor("#888811"),
  68.     plainColor("#888888"), blockColor("#2233CC"), conservedColor("#FF4466");
  70. #define LEFT_JUSTIFY resetiosflags(IOS_BASE::right) << setiosflags(IOS_BASE::left)
  71. #define RIGHT_JUSTIFY resetiosflags(IOS_BASE::left) << setiosflags(IOS_BASE::right)
  73. static string StrToLower(const string& str)
  74. {
  75.     string newStr(str);
  76.     for (int i=0; i<newStr.size(); ++i) newStr[i] = tolower(newStr[i]);
  77.     return newStr;
  78. }
  80. class UnalignedInterval
  81. {
  82. public:
  83.     int alnLocBefore, alnLocAfter, seqLocFrom, seqLocTo;
  84. };
  86. typedef list < UnalignedInterval > IntervalList;
  88. AlignmentDisplay::AlignmentDisplay(const SequenceSet *sequenceSet, const AlignmentSet *alignmentSet) :
  89.     status(CAV_ERROR_DISPLAY)
  90. {
  91.     // start with master row at top, all in lowercase - will be capitalized as things are
  92.     // aligned to it. Also, add an IndexAlnLocToSeqLocRow for the master (only).
  93.     if (!sequenceSet->master) {
  94.         ERR_POST(Error << "Need to know master sequence of SequenceSet for AlignmentDisplay construction");
  95.         return;
  96.     }
  97.     textRows.push_back(new TextRow(StrToLower(sequenceSet->master->sequenceString)));
  99.     // initialize the master index row
  100.     indexAlnLocToSeqLocRows.push_back(
  101.         new IndexAlnLocToSeqLocRow(sequenceSet->master, sequenceSet->master->sequenceString.size()));
  102.     int l;
  103.     for (l=0; l<sequenceSet->master->sequenceString.size(); ++l)
  104.         indexAlnLocToSeqLocRows[0]->SetSeqLocAt(l, l);
  106.     // loop through alignments - one new row for each alignment
  107.     AlignmentSet::AlignmentList::const_iterator a, ae = alignmentSet->alignments.end();
  108.     for (a=alignmentSet->alignments.begin(); a!=ae; ++a) {
  110.         // add index row for each slave that contains only Sequence* for now
  111.         indexAlnLocToSeqLocRows.push_back(new IndexAlnLocToSeqLocRow((*a)->slave));
  113.         // first add blank row that will be filled in with slave
  114.         textRows.push_back(new TextRow(GetWidth()));
  116.         // during this pass through the master, do several things:
  117.         //  - change aligned master residues to uppercase
  118.         //  - fill in aligned slave residues with uppercase
  119.         //  - make a list of all unaligned slave regions, and the current display
  120.         //      coordinates of the aligned slave residues on either side
  121.         int prevAlignedSlaveSeqLoc = -1, prevAlignedSlaveAlnLoc = -1;
  122.         IntervalList intervalList;
  123.         for (int alnLoc=0; alnLoc<=GetWidth(); ++alnLoc) {
  125.             int masterSeqLoc = -1, alignedSlaveSeqLoc = -1;
  127.             if (alnLoc < GetWidth()) {
  128.                 masterSeqLoc = indexAlnLocToSeqLocRows.front()->GetSeqLocAt(alnLoc);
  129.                 if (masterSeqLoc >= 0) {
  130.                     alignedSlaveSeqLoc = (*a)->masterToSlave[masterSeqLoc];
  131.                     if (alignedSlaveSeqLoc >= 0) {
  132.                         textRows.front()->SetCharAt(alnLoc,
  133.                             toupper(textRows.front()->GetCharAt(alnLoc)));
  134.                         textRows.back()->SetCharAt(alnLoc,
  135.                             toupper((*a)->slave->sequenceString[alignedSlaveSeqLoc]));
  136.                     }
  137.                 }
  139.             } else {
  140.                 masterSeqLoc = (*a)->master->sequenceString.size();
  141.                 alignedSlaveSeqLoc = (*a)->slave->sequenceString.size();
  142.             }
  144.             if (alignedSlaveSeqLoc >= 0) {
  145.                 if (alignedSlaveSeqLoc - prevAlignedSlaveSeqLoc > 1) {
  146.                     intervalList.resize(intervalList.size() + 1);
  147.                     intervalList.back().alnLocBefore = prevAlignedSlaveAlnLoc;
  148.                     intervalList.back().alnLocAfter = alnLoc;
  149.                     intervalList.back().seqLocFrom = prevAlignedSlaveSeqLoc + 1;
  150.                     intervalList.back().seqLocTo = alignedSlaveSeqLoc - 1;
  151.                 }
  152.                 prevAlignedSlaveSeqLoc = alignedSlaveSeqLoc;
  153.                 prevAlignedSlaveAlnLoc = alnLoc;
  154.             }
  155.         }
  157.         // now, put in the unaligned regions of the slave. If there is more space
  158.         // than residues, then pad with spaces as necessary. If
  159.         // there isn't enough space, then add gaps to all prior alignments.
  160.         IntervalList::iterator i, ie = intervalList.end();
  161.         int alnLocOffset = 0;   // to track # inserted columns
  162.         for (i=intervalList.begin(); i!=ie; ++i) {
  164.             // compensate for inserted columns in display
  165.             i->alnLocBefore += alnLocOffset;
  166.             i->alnLocAfter += alnLocOffset;
  168.             int
  169.                 displaySpace = i->alnLocAfter - i->alnLocBefore - 1,
  170.                 unalignedLength = i->seqLocTo - i->seqLocFrom + 1,
  171.                 extraSpace = displaySpace - unalignedLength;
  173.             // add gaps to make space
  174.             if (extraSpace < 0) {
  176.                 // where to insert gaps if display space is too small
  177.                 int insertPos;
  178.                 if (i->seqLocFrom == 0) {   // left tail
  179.                     insertPos = 0;
  180.                 } else if (i->seqLocTo == (*a)->slave->sequenceString.size() - 1) { // right tail
  181.                     insertPos = GetWidth();
  182.                 } else {    // between aligned blocks
  183.                     insertPos = i->alnLocAfter - displaySpace / 2;
  184.                 }
  186.                 InsertGaps(-extraSpace, insertPos);
  187.                 alnLocOffset -= extraSpace;
  188.                 extraSpace = 0;
  189.             }
  191.             // fill in unaligned regions with lowercase, right-justifying only for left tail
  192.             for (l=0; l<unalignedLength; ++l) {
  193.                 textRows.back()->SetCharAt(
  194.                     i->alnLocBefore + 1 + ((i->seqLocFrom == 0) ? l + extraSpace : l),
  195.                     tolower((*a)->slave->sequenceString[i->seqLocFrom + l]));
  196.             }
  197.         }
  198.     }
  200.     ERR_POST(Info << "initial alignment display size: " << GetWidth() << " x " << GetNRows());
  202.     // the above algorithm introduces more gaps than are strictly necessary. This
  203.     // will "squeeze" the alignment, deleting gaps wherever possible
  204.     Squeeze();
  206.     // finally, redistribute unaligned residues so that they are split equally between
  207.     // flanking aligned residues
  208.     SplitUnaligned();
  210.     // The Squeeze and SplitUnaligned leave the master index row out of sync, so reindex it
  211.     indexAlnLocToSeqLocRows[0]->ReIndex(*(textRows[0]));
  213.     // find first and last aligned master residues (in alignment coords)
  214.     firstAlnLoc = GetWidth();
  215.     lastAlnLoc = -1;
  216.     for (l=0; l<GetWidth(); ++l)
  217.         if (IsAligned(textRows[0]->GetCharAt(l))) {
  218.             firstAlnLoc = l;
  219.             break;
  220.         }
  221.     for (l=GetWidth()-1; l>=0; --l)
  222.         if (IsAligned(textRows[0]->GetCharAt(l))) {
  223.             lastAlnLoc = l;
  224.             break;
  225.         }
  227.     ERR_POST(Info << "final alignment display size: " << GetWidth() << " x " << GetNRows());
  228.     ERR_POST(Info << "aligned region: " << firstAlnLoc << " through " << lastAlnLoc);
  229.     status = CAV_SUCCESS;
  230. }
  232. AlignmentDisplay::~AlignmentDisplay()
  233. {
  234.     int i;
  235.     for (i=0; i<indexAlnLocToSeqLocRows.size(); ++i) delete indexAlnLocToSeqLocRows[i];
  236.     for (i=0; i<textRows.size(); ++i) delete textRows[i];
  237. }
  239. // shift unaligned residues as far as possible to the left
  240. void AlignmentDisplay::ShiftUnalignedLeft(void)
  241. {
  242.     int gapLoc, resLoc;
  244.     for (int row=0; row<GetNRows(); ++row) {
  246.         gapLoc = 0;
  247.         while (gapLoc < GetWidth()) {
  249.             // find a gap
  250.             while (gapLoc < GetWidth() && !IsGap(textRows[row]->GetCharAt(gapLoc))) ++gapLoc;
  251.             if (gapLoc == GetWidth()) break;
  253.             // find the next unaligned residue
  254.             resLoc = gapLoc + 1;
  255.             while (resLoc < GetWidth() && IsGap(textRows[row]->GetCharAt(resLoc))) ++resLoc;
  256.             if (resLoc == GetWidth()) break;
  257.             if (!IsUnaligned(textRows[row]->GetCharAt(resLoc))) {
  258.                 gapLoc = resLoc + 1;
  259.                 continue;
  260.             }
  262.             // shift unaligned residues over
  263.             while (resLoc < GetWidth() && IsUnaligned(textRows[row]->GetCharAt(resLoc))) {
  264.                 textRows[row]->SetCharAt(gapLoc++, textRows[row]->GetCharAt(resLoc));
  265.                 textRows[row]->SetCharAt(resLoc++, '-');
  266.             }
  267.         }
  268.     }
  269. }
  271. void AlignmentDisplay::Squeeze(void)
  272. {
  273.     // move all unaligned residues as far left as possible - makes it much simpler
  274.     // to identify squeezable regions
  275.     ShiftUnalignedLeft();
  277.     typedef vector < int > SqueezeLocs;
  278.     SqueezeLocs squeezeLocs(GetNRows());
  280.     // find last aligned residue; stop search one after that
  281.     int alnLoc, lastAlignedLoc;
  282.     for (lastAlignedLoc=GetWidth()-2;
  283.          lastAlignedLoc>=0 && !IsAligned(textRows[0]->GetCharAt(lastAlignedLoc));
  284.          --lastAlignedLoc) ;
  285.     ERR_POST(Info << "checking for squeeze up to " << (lastAlignedLoc+1));
  287.     for (alnLoc=0; alnLoc<=lastAlignedLoc+1; ++alnLoc) {
  289.         // check to see whether each row is squeezable at this location
  290.         int row, nGaps, minNGaps = GetWidth();
  291.         for (row=0; row<GetNRows(); ++row) {
  292.             if (!textRows[row]->IsSqueezable(alnLoc, &nGaps, &(squeezeLocs[row]), minNGaps))
  293.                 break;
  294.             if (nGaps < minNGaps) minNGaps = nGaps;
  295.         }
  297.         // if all rows are squeezable, then do the squeeze
  298.         if (row == GetNRows()) {
  299.             ERR_POST(Info << "squeezing " << minNGaps << " gaps at loc " << alnLoc);
  300.             for (row=0; row<GetNRows(); ++row)
  301.                 textRows[row]->DeleteGaps(minNGaps, squeezeLocs[row]);
  302.             lastAlignedLoc -= minNGaps; // account for shift of lastAlignedLoc
  303.         }
  305.         // after checking very first column, skip to first aligned loc to save time
  306.         if (alnLoc == 0)
  307.             while (alnLoc<=lastAlignedLoc && !IsAligned(textRows[0]->GetCharAt(alnLoc))) ++alnLoc;
  308.     }
  309. }
  311. // redistribute unaligned residues so they're split left/right between aligned residues.
  312. // This assumes ShiftUnalignedLeft() has already been called, so that all unaligned
  313. // residues are on the left of the gaps.
  314. void AlignmentDisplay::SplitUnaligned(void)
  315. {
  316.     // alnLocs of various key residues
  317.     int firstAligned, prevAligned, nextAligned, firstUnaligned, lastUnaligned;
  318.     int nGaps, nShift, shiftRes, shiftGap, i;
  320.     for (int row=0; row<GetNRows(); ++row) {
  322.         // find first aligned loc; count gaps up to there
  323.         nGaps = 0;
  324.         for (i=0; i<GetWidth()-2 && !IsAligned(textRows[row]->GetCharAt(i)); ++i)
  325.              if (IsGap(textRows[row]->GetCharAt(i))) ++nGaps;
  326.         firstAligned = i;
  328.         // right-shift left tails
  329.         if (nGaps > 0) {
  330.             for (i=0; i<firstAligned-nGaps; ++i) {
  331.                 textRows[row]->SetCharAt(firstAligned-1-i, textRows[row]->GetCharAt(firstAligned-1-nGaps-i));
  332.                 textRows[row]->SetCharAt(firstAligned-1-nGaps-i, '-');
  333.             }
  334.         }
  336.         prevAligned = firstAligned;
  337.         while (prevAligned < GetWidth()-2) {
  339.             // find aligned res immediately followed by unaligned
  340.             if (!(IsAligned(textRows[row]->GetCharAt(prevAligned)) &&
  341.                   IsUnaligned(textRows[row]->GetCharAt(prevAligned + 1)))) {
  342.                 ++prevAligned;
  343.                 continue;
  344.             }
  345.             firstUnaligned = prevAligned + 1;
  347.             // find last unaligned residue
  348.             for (lastUnaligned = firstUnaligned;
  349.                  lastUnaligned < GetWidth()-1 &&
  350.                     IsUnaligned(textRows[row]->GetCharAt(lastUnaligned + 1));
  351.                  ++lastUnaligned) ;
  353.             // find next aligned after this
  354.             for (nextAligned = lastUnaligned + 1;
  355.                  nextAligned < GetWidth() &&
  356.                     !IsAligned(textRows[row]->GetCharAt(nextAligned));
  357.                  ++nextAligned) ;
  358.             if (nextAligned == GetWidth()) break;
  360.             // shift over the right half of the unaligned stretch
  361.             nGaps = nextAligned - lastUnaligned - 1;
  362.             nShift = (lastUnaligned - firstUnaligned + 1) / 2;
  363.             if (nGaps > 0 && nShift > 0) {
  364.                 shiftRes = lastUnaligned;
  365.                 shiftGap = nextAligned - 1;
  366.                 for (i=0; i<nShift; ++i) {
  367.                     textRows[row]->SetCharAt(shiftGap - i, textRows[row]->GetCharAt(shiftRes - i));
  368.                     textRows[row]->SetCharAt(shiftRes - i, '-');
  369.                 }
  370.             }
  371.             prevAligned = nextAligned;
  372.         }
  373.     }
  374. }
  376. void AlignmentDisplay::InsertGaps(int nGaps, int beforePos)
  377. {
  378.     int i;
  379.     for (i=0; i<indexAlnLocToSeqLocRows.size(); ++i)
  380.         indexAlnLocToSeqLocRows[i]->InsertGaps(nGaps, beforePos);
  381.     for (i=0; i<textRows.size(); ++i)
  382.         textRows[i]->InsertGaps(nGaps, beforePos);
  383. }
  385. char AlignmentDisplay::GetCharAt(int alnLoc, int row) const
  386. {
  387.     if (alnLoc < 0 || alnLoc >= GetWidth() || row < 0 || row >= GetNRows()) {
  388.         ERR_POST(Error << "AlignmentDisplay::GetCharAt() - coordinate out of range");
  389.         return '?';
  390.     }
  392.     return textRows[row]->GetCharAt(alnLoc);
  393. }
  395. class CondensedColumn : public CObject
  396. {
  397. protected:
  398.     const string color;
  399.     vector < int > info;
  401.     CondensedColumn(int nRows, string columnColor) : color(columnColor) { info.resize(nRows, 0); }
  403. public:
  404.     virtual ~CondensedColumn(void) { }
  406.     string GetColor(void) const { return color; }
  407.     virtual int GetDisplayWidth(void) const = 0;
  408.     virtual int GetNResidues(int row) const = 0;
  409.     virtual void DumpRow(CNcbiOstream& os, int row) const = 0;
  410.     virtual void AddRowChar(int row, char ch) = 0;
  411. };
  413. class CondensedColumnAligned : public CondensedColumn
  414. {
  415. public:
  416.     CondensedColumnAligned(int nRows, string color) : CondensedColumn(nRows, color) { }
  417.     int GetDisplayWidth(void) const { return 1; }
  418.     int GetNResidues(int row) const { return 1; }
  419.     void DumpRow(CNcbiOstream& os, int row) const;
  420.     void AddRowChar(int row, char ch);
  421. };
  423. void CondensedColumnAligned::DumpRow(CNcbiOstream& os, int row) const
  424. {
  425.     os << ((char) info[row]);
  426. }
  428. void CondensedColumnAligned::AddRowChar(int row, char ch)
  429. {
  430.     info[row] = (int) ch;
  431. }
  433. class CondensedColumnUnaligned : public CondensedColumn
  434. {
  435. private:
  436.     static const string prefix, postfix;
  437.     int nDigits;
  439. public:
  440.     CondensedColumnUnaligned(int nRows, string color) : CondensedColumn(nRows, color) { nDigits = 1; }
  441.     int GetDisplayWidth(void) const { return prefix.size() + nDigits + postfix.size(); }
  442.     int GetNResidues(int row) const { return info[row]; }
  443.     void DumpRow(CNcbiOstream& os, int row) const;
  444.     void AddRowChar(int row, char ch);
  445. };
  447. const string CondensedColumnUnaligned::prefix = ".[";
  448. const string CondensedColumnUnaligned::postfix = "].";
  450. void CondensedColumnUnaligned::DumpRow(CNcbiOstream& os, int row) const
  451. {
  452.     if (info[row] > 0)
  453.         os << prefix << RIGHT_JUSTIFY << setw(nDigits) << info[row] << postfix;
  454.     else
  455.         os << setw(GetDisplayWidth()) << ' ';
  456. }
  458. void CondensedColumnUnaligned::AddRowChar(int row, char ch)
  459. {
  460.     if (!IsGap(ch)) {
  461.         (info[row])++;
  462.         int digits = ((int) log10((double) info[row])) + 1;
  463.         if (digits > nDigits)
  464.             nDigits = digits;
  465.     }
  466. }
  468. int AlignmentDisplay::DumpCondensed(CNcbiOstream& os, unsigned int options,
  469.     int firstCol, int lastCol, int nColumns, double conservationThreshhold,
  470.     const char *titleHTML, int nFeatures, const AlignmentFeature *alnFeatures) const
  471. {
  472.     bool doHTML = ((options & CAV_HTML) > 0), doHTMLHeader = ((options & CAV_HTML_HEADER) > 0);
  474.     if (firstCol < 0 || lastCol >= GetWidth() || firstCol > lastCol || nColumns < 1) {
  475.         ERR_POST(Error << "AlignmentDisplay::DumpText() - nonsensical display region parameters");
  476.         return CAV_ERROR_BAD_PARAMS;
  477.     }
  479.     // how many rows in the display, including annotations?
  480.     int nDisplayRows = GetNRows();
  481.     if (!alnFeatures) nFeatures = 0;
  482.     if (nFeatures > 0) nDisplayRows += nFeatures;
  484.     // do make look-up location index for each feature
  485.     vector < vector < bool > > annotLocsByIndex;
  486.     int i, featIndex, masterIndex;
  487.     if (nFeatures > 0) {
  488.         annotLocsByIndex.resize(nFeatures);
  489.         for (featIndex=0; featIndex<nFeatures; ++featIndex) {
  490.             annotLocsByIndex[featIndex].resize(indexAlnLocToSeqLocRows[0]->sequence->Length(), false);
  491.             for (i=0; i<alnFeatures[featIndex].nLocations; ++i) {
  492.                 masterIndex = alnFeatures[featIndex].locations[i];
  493.                 if (masterIndex >= 0 && masterIndex < indexAlnLocToSeqLocRows[0]->sequence->Length())
  494.                     annotLocsByIndex[featIndex][masterIndex] = true;
  495.             }
  496.         }
  497.     }
  499.     // condense the display into CondensedColumn list
  500.     deque < CRef < CondensedColumn > > condensedColumns;
  501.     int alnLoc, alnRow, row;
  502.     bool isAlnRow;
  503.     for (alnLoc=firstCol; alnLoc<=lastCol; ++alnLoc) {
  504.         for (alnRow=0; alnRow<GetNRows(); ++alnRow)
  505.             if (!IsAligned(textRows[alnRow]->GetCharAt(alnLoc))) break;
  506.         bool alignedColumn = (alnRow == GetNRows());
  507.         if (alignedColumn) {
  508.             condensedColumns.resize(condensedColumns.size() + 1);
  509.             condensedColumns.back().Reset(
  510.                 new CondensedColumnAligned(nDisplayRows, GetColumnColor(alnLoc, conservationThreshhold)));
  511.         } else {
  512.             if (condensedColumns.size() == 0 ||
  513.                 !(dynamic_cast<CondensedColumnUnaligned*>(condensedColumns.back().GetPointer()))) {
  514.                 condensedColumns.resize(condensedColumns.size() + 1);
  515.                 condensedColumns.back().Reset(
  516.                     new CondensedColumnUnaligned(nDisplayRows, GetColumnColor(alnLoc, conservationThreshhold)));
  517.             }
  518.         }
  519.         for (row=0; row<nDisplayRows; ++row) {
  520.             if (options & CAV_ANNOT_BOTTOM) {
  521.                 alnRow = row;
  522.                 featIndex = row - GetNRows();
  523.             } else {
  524.                 featIndex = row;
  525.                 alnRow = row - nFeatures;
  526.             }
  527.             isAlnRow = (alnRow >= 0 && alnRow < GetNRows());
  528.             if (isAlnRow) {
  529.                 condensedColumns.back()->AddRowChar(row, textRows[alnRow]->GetCharAt(alnLoc));
  530.             } else if (alignedColumn) {
  531.                 // display feature characters only in aligned columns
  532.                 masterIndex = indexAlnLocToSeqLocRows[0]->GetSeqLocAt(alnLoc);
  533.                 if (masterIndex >= 0 && annotLocsByIndex[featIndex][masterIndex])
  534.                     condensedColumns.back()->AddRowChar(row, alnFeatures[featIndex].featChar);
  535.                 else
  536.                     condensedColumns.back()->AddRowChar(row, ' ');
  537.             }
  538.         }
  539.     }
  541.     // set up the titles and uids, figure out how much space any seqLoc string will take
  542.     vector < string > titles(nDisplayRows), uids(doHTML ? GetNRows() : 0);
  543.     int maxTitleLength = 0, maxSeqLocStrLength = 0, leftMargin, decimalLength;
  544.     for (row=0; row<nDisplayRows; ++row) {
  545.         if (options & CAV_ANNOT_BOTTOM) {
  546.             alnRow = row;
  547.             featIndex = row - GetNRows();
  548.         } else {
  549.             featIndex = row;
  550.             alnRow = row - nFeatures;
  551.         }
  552.         isAlnRow = (alnRow >= 0 && alnRow < GetNRows());
  554.         // title
  555.         titles[row] = isAlnRow ?
  556.             indexAlnLocToSeqLocRows[alnRow]->sequence->GetTitle() :
  557.             string(alnFeatures[featIndex].shortName);
  558.         if (titles[row].size() > maxTitleLength) maxTitleLength = titles[row].size();
  559.         if (isAlnRow) {
  560.             decimalLength = ((int) log10((double)
  561.                 indexAlnLocToSeqLocRows[alnRow]->sequence->sequenceString.size())) + 1;
  562.             if (decimalLength > maxSeqLocStrLength) maxSeqLocStrLength = decimalLength;
  563.         }
  565.         // uid for link to entrez
  566.         if (doHTML && isAlnRow) {
  567.             if (indexAlnLocToSeqLocRows[alnRow]->sequence->pdbID.size() > 0) {
  568.                 if (indexAlnLocToSeqLocRows[alnRow]->sequence->pdbID != "query" &&
  569.                     indexAlnLocToSeqLocRows[alnRow]->sequence->pdbID != "consensus") {
  570.                     uids[alnRow] = indexAlnLocToSeqLocRows[alnRow]->sequence->pdbID;
  571.                     if (indexAlnLocToSeqLocRows[alnRow]->sequence->pdbChain != ' ')
  572.                         uids[alnRow] += (char) indexAlnLocToSeqLocRows[alnRow]->sequence->pdbChain;
  573.                 }
  574.             } else if (indexAlnLocToSeqLocRows[alnRow]->sequence->gi != Sequence::NOT_SET) {
  575.                 CNcbiOstrstream uidoss;
  576.                 uidoss << indexAlnLocToSeqLocRows[alnRow]->sequence->gi << '';
  577.                 uids[alnRow] = uidoss.str();
  578.                 delete uidoss.str();
  579.             } else if (indexAlnLocToSeqLocRows[alnRow]->sequence->accession.size() > 0) {
  580.                 uids[alnRow] = indexAlnLocToSeqLocRows[alnRow]->sequence->accession;
  581.             }
  582.         }
  583.     }
  584.     leftMargin = maxTitleLength + maxSeqLocStrLength + 2;
  586.     // need to keep track of first, last seqLocs for each row in each paragraph;
  587.     // find seqLoc of first residue >= firstCol
  588.     vector < int > lastShownSeqLocs(GetNRows());
  589.     for (alnRow=0; alnRow<GetNRows(); ++alnRow) {
  590.         lastShownSeqLocs[alnRow] = -1;
  591.         for (alnLoc=0; alnLoc<firstCol; ++alnLoc)
  592.             if (!IsGap(textRows[alnRow]->GetCharAt(alnLoc))) lastShownSeqLocs[alnRow]++;
  593.     }
  595.     // header
  596.     if (doHTML && doHTMLHeader)
  597.         os << "<HTML><TITLE>" << (titleHTML ? titleHTML : "CDDAlignView HTML Display") <<
  598.             "</TITLE><BODY BGCOLOR=" << bgColor << ">n";;
  600.     // split alignment up into "paragraphs", each with <= nColumns
  601.     if (doHTML) os << "<TABLE>n";
  602.     int paragraphStart, nParags = 0, nCondensedColumns;
  603.     ERR_POST(Info << "paragraph width: " << nColumns);
  604.     for (paragraphStart=0;
  605.          paragraphStart<condensedColumns.size();
  606.          paragraphStart+=nCondensedColumns, ++nParags) {
  608.         // figure out how many condensed columns will fit in this paragraph
  609.         int displayWidth = condensedColumns[paragraphStart]->GetDisplayWidth();
  610.         nCondensedColumns = 1;
  611.         while (paragraphStart+nCondensedColumns < condensedColumns.size()) {
  612.             int columnWidth = condensedColumns[paragraphStart+nCondensedColumns]->GetDisplayWidth();
  613.             if (displayWidth + columnWidth <= nColumns) {
  614.                 displayWidth += columnWidth;
  615.                 ++nCondensedColumns;
  616.             } else
  617.                 break;
  618.         }
  620.         // start table row
  621.         if (doHTML)
  622.             os << "<tr><td bgcolor=" << blockBGColors[nParags % nBlockColors] << "><pre>nn";
  623.         else
  624.             if (paragraphStart > 0) os << 'n';
  626.         // output each alignment row
  627.         for (row=0; row<nDisplayRows; ++row) {
  628.             if (options & CAV_ANNOT_BOTTOM) {
  629.                 alnRow = row;
  630.                 featIndex = row - GetNRows();
  631.             } else {
  632.                 featIndex = row;
  633.                 alnRow = row - nFeatures;
  634.             }
  635.             isAlnRow = (alnRow >= 0 && alnRow < GetNRows());
  637.             // title
  638.             if (isAlnRow && doHTML && uids[alnRow].size() > 0) {
  639.                 os << "<a href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi"
  640.                     << "?cmd=Search&doptcmdl=GenPept&db=Protein&term="
  641.                     << uids[alnRow] << "" onMouseOut="window.status=''"n"
  642.                     << "onMouseOver="window.status='"
  643.                     << ((indexAlnLocToSeqLocRows[alnRow]->sequence->description.size() > 0) ?
  644.                             indexAlnLocToSeqLocRows[alnRow]->sequence->description : titles[row])
  645.                     << "';return true">"
  646.                     << setw(0) << titles[row] << "</a>";
  647.             } else {
  648.                 if (doHTML && !isAlnRow) os << "<font color=" << featColor << '>';
  649.                 os << setw(0) << titles[row];
  650.             }
  651.             os << setw(maxTitleLength+1-titles[row].size()) << ' ';
  653.             if (isAlnRow) {
  654.                 // count displayed residues
  655.                 int nDisplayedResidues = 0;
  656.                 for (i=0; i<nCondensedColumns; ++i)
  657.                     nDisplayedResidues += condensedColumns[paragraphStart+i]->GetNResidues(row);
  659.                 // left start pos (output 1-numbered for humans...)
  660.                 if (doHTML) os << "<font color=" << numColor << '>';
  661.                 if (nDisplayedResidues > 0)
  662.                     os << RIGHT_JUSTIFY << setw(maxSeqLocStrLength) << (lastShownSeqLocs[alnRow]+2) << ' ';
  663.                 else
  664.                     os << RIGHT_JUSTIFY << setw(maxSeqLocStrLength) << ' ' << ' ';
  666.                 // dump sequence, applying color changes only when necessary
  667.                 if (doHTML) {
  668.                     string prevColor;
  669.                     for (i=0; i<nCondensedColumns; ++i) {
  670.                         string color = condensedColumns[paragraphStart+i]->GetColor();
  671.                         if (color != prevColor) {
  672.                             os << "</font><font color=" << color << '>';
  673.                             prevColor = color;
  674.                         }
  675.                         condensedColumns[paragraphStart+i]->DumpRow(os, row);
  676.                     }
  677.                     os << "</font>";
  678.                 } else {
  679.                     for (i=0; i<nCondensedColumns; ++i) {
  680.                         condensedColumns[paragraphStart+i]->DumpRow(os, row);
  681.                     }
  682.                 }
  684.                 // right end pos
  685.                 if (nDisplayedResidues > 0) {
  686.                     os << ' ';
  687.                     if (doHTML) os << "<font color=" << numColor << '>';
  688.                     os << LEFT_JUSTIFY << setw(0) << (lastShownSeqLocs[alnRow]+nDisplayedResidues+1);
  689.                     if (doHTML) os << "</font>";
  690.                 }
  691.                 os << 'n';
  693.                 // setup to begin next parag
  694.                 lastShownSeqLocs[alnRow] += nDisplayedResidues;
  695.             }
  697.             // print alignment annotation characters
  698.             else {
  699.                 // skip number
  700.                 os << RIGHT_JUSTIFY << setw(maxSeqLocStrLength) << ' ' << ' ';
  702.                 // do characters
  703.                 for (i=0; i<nCondensedColumns; ++i)
  704.                     condensedColumns[paragraphStart+i]->DumpRow(os, row);
  706.                 os << (doHTML ? "</font>n" : "n");
  707.             }
  708.         }
  710.         // end table row
  711.         if (doHTML) os << "</pre></td></tr>n";
  712.     }
  714.     if (doHTML) os << "</TABLE>n";
  716.     // add feature legend
  717.     if (nFeatures > 0) {
  718.         os << (doHTML ? "<BR>n" : "n");
  719.         for (featIndex=0; featIndex<nFeatures; ++featIndex)
  720.             if (alnFeatures[featIndex].description)
  721.                 os << alnFeatures[featIndex].shortName << ": " << alnFeatures[featIndex].description
  722.                    << (doHTML ? "<BR>n" : "n");
  723.     }
  725.     if (doHTML && doHTMLHeader) os << "</BODY></HTML>n";
  727.     // additional sanity check on seqloc markers
  728.     if (firstCol == 0 && lastCol == GetWidth()-1) {
  729.         for (alnRow=0; alnRow<GetNRows(); ++alnRow) {
  730.             if (lastShownSeqLocs[alnRow] !=
  731.                     indexAlnLocToSeqLocRows[alnRow]->sequence->sequenceString.size()-1) {
  732.                 ERR_POST(Error << "full display - seqloc markers don't add up");
  733.                 break;
  734.             }
  735.         }
  736.         if (alnRow == GetNRows())
  737.             ERR_POST(Info << "full display - seqloc markers add up correctly");
  738.     }
  740.     return CAV_SUCCESS;
  741. }
  743. int AlignmentDisplay::DumpText(CNcbiOstream& os, unsigned int options,
  744.     int firstCol, int lastCol, int nColumns, double conservationThreshhold,
  745.     const char *titleHTML, int nFeatures, const AlignmentFeature *alnFeatures) const
  746. {
  747.     bool doHTML = ((options & CAV_HTML) > 0), doHTMLHeader = ((options & CAV_HTML_HEADER) > 0);
  749.     if (firstCol < 0 || lastCol >= GetWidth() || firstCol > lastCol || nColumns < 1) {
  750.         ERR_POST(Error << "AlignmentDisplay::DumpText() - nonsensical display region parameters");
  751.         return CAV_ERROR_BAD_PARAMS;
  752.     }
  754.     // how many rows in the display, including annotations?
  755.     int nDisplayRows = GetNRows();
  756.     if (!alnFeatures) nFeatures = 0;
  757.     if (nFeatures > 0) nDisplayRows += nFeatures;
  759.     // set up the titles and uids, figure out how much space any seqLoc string will take
  760.     vector < string > titles(nDisplayRows), uids(doHTML ? GetNRows() : 0);
  761.     int row, featIndex, alnRow, maxTitleLength = 0, maxSeqLocStrLength = 0, leftMargin, decimalLength;
  762.     bool isAlnRow;
  763.     for (row=0; row<nDisplayRows; ++row) {
  764.         if (options & CAV_ANNOT_BOTTOM) {
  765.             alnRow = row;
  766.             featIndex = row - GetNRows();
  767.         } else {
  768.             featIndex = row;
  769.             alnRow = row - nFeatures;
  770.         }
  771.         isAlnRow = (alnRow >= 0 && alnRow < GetNRows());
  773.         // title
  774.         titles[row] = isAlnRow ?
  775.             indexAlnLocToSeqLocRows[alnRow]->sequence->GetTitle() :
  776.             string(alnFeatures[featIndex].shortName);
  777.         if (titles[row].size() > maxTitleLength) maxTitleLength = titles[row].size();
  778.         if (isAlnRow) {
  779.             decimalLength = ((int) log10((double)
  780.                 indexAlnLocToSeqLocRows[alnRow]->sequence->sequenceString.size())) + 1;
  781.             if (decimalLength > maxSeqLocStrLength) maxSeqLocStrLength = decimalLength;
  782.         }
  784.         // uid for link to entrez
  785.         if (doHTML && isAlnRow) {
  786.             if (indexAlnLocToSeqLocRows[alnRow]->sequence->pdbID.size() > 0) {
  787.                 if (indexAlnLocToSeqLocRows[alnRow]->sequence->pdbID != "query" &&
  788.                     indexAlnLocToSeqLocRows[alnRow]->sequence->pdbID != "consensus") {
  789.                     uids[alnRow] = indexAlnLocToSeqLocRows[alnRow]->sequence->pdbID;
  790.                     if (indexAlnLocToSeqLocRows[alnRow]->sequence->pdbChain != ' ')
  791.                         uids[alnRow] += (char) indexAlnLocToSeqLocRows[alnRow]->sequence->pdbChain;
  792.                 }
  793.             } else if (indexAlnLocToSeqLocRows[alnRow]->sequence->gi != Sequence::NOT_SET) {
  794.                 CNcbiOstrstream uidoss;
  795.                 uidoss << indexAlnLocToSeqLocRows[alnRow]->sequence->gi << '';
  796.                 uids[alnRow] = uidoss.str();
  797.                 delete uidoss.str();
  798.             } else if (indexAlnLocToSeqLocRows[alnRow]->sequence->accession.size() > 0) {
  799.                 uids[alnRow] = indexAlnLocToSeqLocRows[alnRow]->sequence->accession;
  800.             }
  801.         }
  802.     }
  803.     leftMargin = maxTitleLength + maxSeqLocStrLength + 2;
  805.     // need to keep track of first, last seqLocs for each row in each paragraph;
  806.     // find seqLoc of first residue >= firstCol
  807.     vector < int > lastShownSeqLocs(GetNRows());
  808.     int alnLoc;
  809.     for (alnRow=0; alnRow<GetNRows(); ++alnRow) {
  810.         lastShownSeqLocs[alnRow] = -1;
  811.         for (alnLoc=0; alnLoc<firstCol; ++alnLoc)
  812.             if (!IsGap(textRows[alnRow]->GetCharAt(alnLoc))) lastShownSeqLocs[alnRow]++;
  813.     }
  815.     // header
  816.     if (doHTML && doHTMLHeader)
  817.         os << "<HTML><TITLE>" << (titleHTML ? titleHTML : "CDDAlignView HTML Display") <<
  818.             "</TITLE><BODY BGCOLOR=" << bgColor << ">n";;
  820.     // do make look-up location index for each feature
  821.     vector < vector < bool > > annotLocsByIndex;
  822.     int i, masterIndex;
  823.     if (nFeatures > 0) {
  824.         annotLocsByIndex.resize(nFeatures);
  825.         for (featIndex=0; featIndex<nFeatures; ++featIndex) {
  826.             annotLocsByIndex[featIndex].resize(indexAlnLocToSeqLocRows[0]->sequence->Length(), false);
  827.             for (i=0; i<alnFeatures[featIndex].nLocations; ++i) {
  828.                 masterIndex = alnFeatures[featIndex].locations[i];
  829.                 if (masterIndex >= 0 && masterIndex < indexAlnLocToSeqLocRows[0]->sequence->Length())
  830.                     annotLocsByIndex[featIndex][masterIndex] = true;
  831.             }
  832.         }
  833.     }
  835.     // split alignment up into "paragraphs", each with nColumns
  836.     if (doHTML) os << "<TABLE>n";
  837.     int paragraphStart, nParags = 0;
  838.     for (paragraphStart=0; (firstCol+paragraphStart)<=lastCol; paragraphStart+=nColumns, ++nParags) {
  840.         // start table row
  841.         if (doHTML)
  842.             os << "<tr><td bgcolor=" << blockBGColors[nParags % nBlockColors] << "><pre>n";
  843.         else
  844.             if (paragraphStart > 0) os << 'n';
  846.         // do ruler
  847.         int nMarkers = 0, width;
  848.         if (doHTML) os << "<font color=" << rulerColor << '>';
  849.         for (i=0; i<nColumns && (firstCol+paragraphStart+i)<=lastCol; ++i) {
  850.             if ((paragraphStart+i+1)%10 == 0) {
  851.                 if (nMarkers == 0)
  852.                     width = leftMargin + i + 1;
  853.                 else
  854.                     width = 10;
  855.                 os << RIGHT_JUSTIFY << setw(width) << (paragraphStart+i+1);
  856.                 ++nMarkers;
  857.             }
  858.         }
  859.         if (doHTML) os << "</font>";
  860.         os << 'n';
  861.         if (doHTML) os << "<font color=" << rulerColor << '>';
  862.         for (i=0; i<leftMargin; ++i) os << ' ';
  863.         for (i=0; i<nColumns && (firstCol+paragraphStart+i)<=lastCol; ++i) {
  864.             if ((paragraphStart+i+1)%10 == 0)
  865.                 os << '|';
  866.             else if ((paragraphStart+i+1)%5 == 0)
  867.                 os << '*';
  868.             else
  869.                 os << '.';
  870.         }
  871.         if (doHTML) os << "</font>";
  872.         os << 'n';
  874.         // get column colors
  875.         vector < string > columnColors;
  876.         if (doHTML)
  877.             for (i=0; i<nColumns && (firstCol+paragraphStart+i)<=lastCol; ++i)
  878.                 columnColors.resize(columnColors.size()+1,
  879.                     GetColumnColor(firstCol+paragraphStart+i, conservationThreshhold));
  881.         // output each alignment row
  882.         int nDisplayedResidues;
  883.         for (row=0; row<nDisplayRows; ++row) {
  884.             if (options & CAV_ANNOT_BOTTOM) {
  885.                 alnRow = row;
  886.                 featIndex = row - GetNRows();
  887.             } else {
  888.                 featIndex = row;
  889.                 alnRow = row - nFeatures;
  890.             }
  891.             isAlnRow = (alnRow >= 0 && alnRow < GetNRows());
  893.             // actual sequence characters; count how many non-gaps in each row
  894.             nDisplayedResidues = 0;
  895.             string rowChars;
  896.             if (isAlnRow) {
  897.                 for (i=0; i<nColumns && (firstCol+paragraphStart+i)<=lastCol; ++i) {
  898.                     char ch = textRows[alnRow]->GetCharAt(firstCol+paragraphStart+i);
  899.                     rowChars += ch;
  900.                     if (!IsGap(ch)) ++nDisplayedResidues;
  901.                 }
  902.             }
  904.             // title
  905.             if (isAlnRow && doHTML && uids[alnRow].size() > 0) {
  906.                 os << "<a href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi"
  907.                     << "?cmd=Search&doptcmdl=GenPept&db=Protein&term="
  908.                     << uids[alnRow] << "" onMouseOut="window.status=''"n"
  909.                     << "onMouseOver="window.status='"
  910.                     << ((indexAlnLocToSeqLocRows[alnRow]->sequence->description.size() > 0) ?
  911.                             indexAlnLocToSeqLocRows[alnRow]->sequence->description : titles[row])
  912.                     << "';return true">"
  913.                     << setw(0) << titles[row] << "</a>";
  914.             } else {
  915.                 if (doHTML && !isAlnRow) os << "<font color=" << featColor << '>';
  916.                 os << setw(0) << titles[row];
  917.             }
  918.             os << setw(maxTitleLength+1-titles[row].size()) << ' ';
  920.             if (isAlnRow) {
  921.                 // left start pos (output 1-numbered for humans...)
  922.                 if (doHTML) os << "<font color=" << numColor << '>';
  923.                 if (nDisplayedResidues > 0)
  924.                     os << RIGHT_JUSTIFY << setw(maxSeqLocStrLength) << (lastShownSeqLocs[alnRow]+2) << ' ';
  925.                 else
  926.                     os << RIGHT_JUSTIFY << setw(maxSeqLocStrLength) << ' ' << ' ';
  928.                 // dump sequence, applying color changes only when necessary
  929.                 if (doHTML) {
  930.                     string prevColor;
  931.                     for (i=0; i<rowChars.size(); ++i) {
  932.                         if (columnColors[i] != prevColor) {
  933.                             os << "</font><font color=" << columnColors[i] << '>';
  934.                             prevColor = columnColors[i];
  935.                         }
  936.                         os << rowChars[i];
  937.                     }
  938.                     os << "</font>";
  939.                 } else
  940.                     os << rowChars;
  942.                 // right end pos
  943.                 if (nDisplayedResidues > 0) {
  944.                     os << ' ';
  945.                     if (doHTML) os << "<font color=" << numColor << '>';
  946.                     os << LEFT_JUSTIFY << setw(0) << (lastShownSeqLocs[alnRow]+nDisplayedResidues+1);
  947.                     if (doHTML) os << "</font>";
  948.                 }
  949.                 os << 'n';
  951.                 // setup to begin next parag
  952.                 lastShownSeqLocs[alnRow] += nDisplayedResidues;
  953.             }
  955.             // print alignment annotation characters
  956.             else {
  957.                 // skip number
  958.                 os << RIGHT_JUSTIFY << setw(maxSeqLocStrLength) << ' ' << ' ';
  960.                 // do characters, but only allow annot where master residue is aligned to something
  961.                 for (i=0; i<nColumns && (firstCol+paragraphStart+i)<=lastCol; ++i) {
  962.                     if (IsAligned(textRows[0]->GetCharAt(firstCol+paragraphStart+i))) {
  963.                         masterIndex = indexAlnLocToSeqLocRows[0]->GetSeqLocAt(firstCol+paragraphStart+i);
  964.                         os << ((masterIndex >= 0 && annotLocsByIndex[featIndex][masterIndex])
  965.                             ? alnFeatures[featIndex].featChar : ' ');
  966.                     } else
  967.                         os << ' ';
  968.                 }
  970.                 os << (doHTML ? "</font>n" : "n");
  971.             }
  972.         }
  974.         // end table row
  975.         if (doHTML) os << "</pre></td></tr>n";
  976.     }
  978.     if (doHTML) os << "</TABLE>n";
  980.     // add feature legend
  981.     if (nFeatures > 0) {
  982.         os << (doHTML ? "<BR>n" : "n");
  983.         for (featIndex=0; featIndex<nFeatures; ++featIndex)
  984.             if (alnFeatures[featIndex].description)
  985.                 os << alnFeatures[featIndex].shortName << ": " << alnFeatures[featIndex].description
  986.                    << (doHTML ? "<BR>n" : "n");
  987.     }
  989.     if (doHTML && doHTMLHeader) os << "</BODY></HTML>n";
  991.     // additional sanity check on seqloc markers
  992.     if (firstCol == 0 && lastCol == GetWidth()-1) {
  993.         for (alnRow=0; alnRow<GetNRows(); ++alnRow) {
  994.             if (lastShownSeqLocs[alnRow] !=
  995.                     indexAlnLocToSeqLocRows[alnRow]->sequence->sequenceString.size()-1) {
  996.                 ERR_POST(Error << "full display - seqloc markers don't add up");
  997.                 break;
  998.             }
  999.         }
  1000.         if (alnRow == GetNRows())
  1001.             ERR_POST(Info << "full display - seqloc markers add up correctly");
  1002.     }
  1004.     return CAV_SUCCESS;
  1005. }
  1007. int AlignmentDisplay::DumpFASTA(int firstCol, int lastCol, int nColumns,
  1008.     bool doLowercase, CNcbiOstream& os) const
  1009. {
  1010.     if (firstCol < 0 || lastCol >= GetWidth() || firstCol > lastCol || nColumns < 1) {
  1011.         ERR_POST(Error << "AlignmentDisplay::DumpFASTA() - nonsensical display region parameters");
  1012.         return CAV_ERROR_BAD_PARAMS;
  1013.     }
  1015.     // output each alignment row
  1016.     for (int row=0; row<GetNRows(); ++row) {
  1018.         // create title line
  1019.         os << '>';
  1020.         const Sequence *seq = indexAlnLocToSeqLocRows[row]->sequence;
  1021.         bool prevID = false;
  1022.         if (seq->gi != Sequence::NOT_SET) {
  1023.             os << "gi|" << seq->gi;
  1024.             prevID = true;
  1025.         }
  1026.         if (seq->pdbID.size() > 0) {
  1027.             if (prevID) os << '|';
  1028.             if (seq->pdbID == "query" || seq->pdbID == "consensus") {
  1029.                 os << "lcl|" << seq->pdbID;
  1030.             } else {
  1031.                 os << "pdb|" << seq->pdbID;
  1032.                 if (seq->pdbChain != ' ')
  1033.                     os << '|' << (char) seq->pdbChain << " Chain "
  1034.                        << (char) seq->pdbChain << ',';
  1035.             }
  1036.             prevID = true;
  1037.         }
  1038.         if (seq->accession.size() > 0) {
  1039.             if (prevID) os << '|';
  1040.             os << seq->accession;
  1041.             prevID = true;
  1042.         }
  1043.         if (seq->description.size() > 0)
  1044.             os << ' ' << seq->description;
  1045.         os << 'n';
  1047.         // split alignment up into "paragraphs", each with nColumns
  1048.         int paragraphStart, nParags = 0, i;
  1049.         for (paragraphStart=0; (firstCol+paragraphStart)<=lastCol; paragraphStart+=nColumns, ++nParags) {
  1050.             for (i=0; i<nColumns && (firstCol+paragraphStart+i)<=lastCol; ++i) {
  1051.                 char ch = textRows[row]->GetCharAt(firstCol+paragraphStart+i);
  1052.                 if (!doLowercase) ch = toupper(ch);
  1053.                 os << ch;
  1054.             }
  1055.             os << 'n';
  1056.         }
  1057.     }
  1059.     return CAV_SUCCESS;
  1060. }
  1062. const string& AlignmentDisplay::GetColumnColor(int alnLoc, double conservationThreshhold) const
  1063. {
  1064.     // standard probabilities (calculated by BLAST using BLOSUM62 - see conservation_colorer.cpp in Cn3D++)
  1065.     typedef map < char , double > Char2Double;
  1066.     static Char2Double StandardProbabilities;
  1067.     if (StandardProbabilities.size() == 0) {
  1068.         StandardProbabilities['A'] = 0.07805;
  1069.         StandardProbabilities['C'] = 0.01925;
  1070.         StandardProbabilities['D'] = 0.05364;
  1071.         StandardProbabilities['E'] = 0.06295;
  1072.         StandardProbabilities['F'] = 0.03856;
  1073.         StandardProbabilities['G'] = 0.07377;
  1074.         StandardProbabilities['H'] = 0.02199;
  1075.         StandardProbabilities['I'] = 0.05142;
  1076.         StandardProbabilities['K'] = 0.05744;
  1077.         StandardProbabilities['L'] = 0.09019;
  1078.         StandardProbabilities['M'] = 0.02243;
  1079.         StandardProbabilities['N'] = 0.04487;
  1080.         StandardProbabilities['P'] = 0.05203;
  1081.         StandardProbabilities['Q'] = 0.04264;
  1082.         StandardProbabilities['R'] = 0.05129;
  1083.         StandardProbabilities['S'] = 0.07120;
  1084.         StandardProbabilities['T'] = 0.05841;
  1085.         StandardProbabilities['V'] = 0.06441;
  1086.         StandardProbabilities['W'] = 0.01330;
  1087.         StandardProbabilities['X'] = 0;
  1088.         StandardProbabilities['Y'] = 0.03216;
  1089.     }
  1091.     // if this column isn't completely aligned, use plain color
  1092.     int row;
  1093.     for (row=0; row<GetNRows(); ++row)
  1094.         if (!IsAligned(textRows[row]->GetCharAt(alnLoc))) return plainColor;
  1096.     // create profile (residue frequencies) for this column
  1097.     Char2Double profile;
  1098.     Char2Double::iterator p, pe;
  1099.     for (row=0; row<GetNRows(); ++row) {
  1100.         char ch = toupper(textRows[row]->GetCharAt(alnLoc));
  1101.         switch (ch) {
  1102.             case 'A': case 'R': case 'N': case 'D': case 'C':
  1103.             case 'Q': case 'E': case 'G': case 'H': case 'I':
  1104.             case 'L': case 'K': case 'M': case 'F': case 'P':
  1105.             case 'S': case 'T': case 'W': case 'Y': case 'V':
  1106.                 break;
  1107.             default:
  1108.                 ch = 'X'; // make all but natural aa's just 'X'
  1109.         }
  1110.         if ((p=profile.find(ch)) != profile.end())
  1111.             p->second += 1.0/GetNRows();
  1112.         else
  1113.             profile[ch] = 1.0/GetNRows();
  1114.     }
  1116.     // check for identity...
  1117.     if (conservationThreshhold == SHOW_IDENTITY) {
  1118.         if (profile.size() == 1 && profile.begin()->first != 'X')
  1119.             return conservedColor;
  1120.         else
  1121.             return blockColor;
  1122.     }
  1124.     // ... or calculate information content for this column (calculated in bits -> logs of base 2)
  1125.     double information = 0.0;
  1126.     pe = profile.end();
  1127.     for (p=profile.begin(); p!=pe; ++p) {
  1128.         static const double ln2 = log(2.0), threshhold = 0.0001;
  1129.         double expFreq = StandardProbabilities[p->first];
  1130.         if (expFreq > threshhold) {
  1131.             float freqRatio = p->second / expFreq;
  1132.             if (freqRatio > threshhold)
  1133.                 information += p->second * log(freqRatio) / ln2;
  1134.         }
  1135.     }
  1137.     // if conserved, use conservation color
  1138.     if (information > conservationThreshhold) return conservedColor;
  1140.     // if is unconserved, use block color to show post-IBM block locations
  1141.     return blockColor;
  1142. }
  1145. ///// Row class methods /////
  1147. void TextRow::InsertGaps(int nGaps, int beforePos)
  1148. {
  1149.     if (beforePos < 0 || beforePos > Length()) {
  1150.         ERR_POST(Error << "TextRow::InsertGaps() - beforePos out of range");
  1151.         return;
  1152.     }
  1154.     chars.insert(beforePos, nGaps, '-');
  1155. }
  1157. void TextRow::DeleteGaps(int nGaps, int startPos)
  1158. {
  1159.     if (startPos < 0 || startPos+nGaps-1 > Length()) {
  1160.         ERR_POST(Error << "TextRow::DeleteGaps() - startPos out of range");
  1161.         return;
  1162.     }
  1164.     // check to make sure they're all gaps
  1165.     for (int i=0; i<nGaps; ++i) {
  1166.         if (!IsGap(chars[startPos + i])) {
  1167.             ERR_POST(Error << "TextRow::DeleteGaps() - trying to delete non-gap");
  1168.             return;
  1169.         }
  1170.     }
  1172.     chars.erase(startPos, nGaps);
  1173. }
  1175. // find out how many gaps (up to maxGaps) are present from alnLoc to the next aligned
  1176. // residue to the right; set startPos to the first gap, nGaps to # gaps starting at startPos
  1177. bool TextRow::IsSqueezable(int alnLoc, int *nGaps, int *startPos, int maxGaps) const
  1178. {
  1179.     if (alnLoc < 0 || alnLoc >= chars.size()) {
  1180.         ERR_POST(Error << "TextRow::IsSqueezable() - alnLoc out of range");
  1181.         return false;
  1182.     }
  1184.     // skip unaligned residues
  1185.     while (alnLoc < chars.size() && IsUnaligned(chars[alnLoc])) ++alnLoc;
  1186.     if (alnLoc == chars.size() || IsAligned(chars[alnLoc])) return false;
  1188.     // count gaps
  1189.     *startPos = alnLoc;
  1190.     for (*nGaps=1, ++alnLoc; alnLoc < chars.size() && IsGap(chars[alnLoc]); (*nGaps)++, ++alnLoc)
  1191.         if (*nGaps == maxGaps) break;
  1192.     return true;
  1193. }
  1196. IndexAlnLocToSeqLocRow::IndexAlnLocToSeqLocRow(const Sequence *seq, int length) :
  1197.     sequence(seq)
  1198. {
  1199.     if (length > 0) seqLocs.resize(length, -1);
  1200. }
  1202. void IndexAlnLocToSeqLocRow::InsertGaps(int nGaps, int beforePos)
  1203. {
  1204.     if (nGaps <= 0 || seqLocs.size() == 0) return;
  1206.     if (beforePos < 0 || beforePos > Length()) {
  1207.         ERR_POST(Error << "IndexAlnLocToSeqLocRow::InsertGaps() - beforePos out of range");
  1208.         return;
  1209.     }
  1211.     IntVec::iterator s = seqLocs.begin();
  1212.     for (int i=0; i<beforePos; ++i) ++s;
  1213.     seqLocs.insert(s, nGaps, -1);
  1214. }
  1216. void IndexAlnLocToSeqLocRow::ReIndex(const TextRow& textRow)
  1217. {
  1218.     seqLocs.resize(textRow.Length());
  1219.     int seqLoc = 0;
  1220.     for (int i=0; i<textRow.Length(); ++i) {
  1221.         if (IsGap(textRow.GetCharAt(i)))
  1222.             seqLocs[i] = -1;
  1223.         else
  1224.             seqLocs[i] = seqLoc++;
  1225.     }
  1226.     if (seqLoc != sequence->sequenceString.size())
  1227.         ERR_POST(Error << "IndexAlnLocToSeqLocRow::ReIndex() - wrong sequence length");
  1228. }
  1230. END_NCBI_SCOPE
  1233. /*
  1234. * ---------------------------------------------------------------------------
  1235. * $Log: cav_alndisplay.cpp,v $
  1236. * Revision 1000.3  2004/06/01 19:41:15  gouriano
  1237. * PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.5
  1238. *
  1239. * Revision 1.5  2004/05/21 21:42:51  gorelenk
  1240. * Added PCH ncbi_pch.hpp
  1241. *
  1242. * Revision 1.4  2004/03/15 18:51:27  thiessen
  1243. * prefer prefix vs. postfix ++/-- operators
  1244. *
  1245. * Revision 1.3  2003/11/15 13:16:04  thiessen
  1246. * fix uid link urls
  1247. *
  1248. * Revision 1.2  2003/06/02 16:06:41  dicuccio
  1249. * Rearranged src/objects/ subtree.  This includes the following shifts:
  1250. *     - src/objects/asn2asn --> arc/app/asn2asn
  1251. *     - src/objects/testmedline --> src/objects/ncbimime/test
  1252. *     - src/objects/objmgr --> src/objmgr
  1253. *     - src/objects/util --> src/objmgr/util
  1254. *     - src/objects/alnmgr --> src/objtools/alnmgr
  1255. *     - src/objects/flat --> src/objtools/flat
  1256. *     - src/objects/validator --> src/objtools/validator
  1257. *     - src/objects/cddalignview --> src/objtools/cddalignview
  1258. * In addition, libseq now includes six of the objects/seq... libs, and libmmdb
  1259. * replaces the three libmmdb? libs.
  1260. *
  1261. * Revision 1.1  2003/03/19 19:04:12  thiessen
  1262. * move again
  1263. *
  1264. * Revision 1.2  2003/03/19 13:25:09  thiessen
  1265. * fix log ambiguity
  1266. *
  1267. * Revision 1.1  2003/03/19 05:33:43  thiessen
  1268. * move to src/app/cddalignview
  1269. *
  1270. * Revision 1.23  2003/02/03 17:52:03  thiessen
  1271. * move CVS Log to end of file
  1272. *
  1273. * Revision 1.22  2003/01/21 18:01:07  thiessen
  1274. * add condensed alignment display
  1275. *
  1276. * Revision 1.21  2002/12/09 13:31:04  thiessen
  1277. * use new query.fcgi for genpept links
  1278. *
  1279. * Revision 1.20  2002/11/08 19:38:11  thiessen
  1280. * add option for lowercase unaligned in FASTA
  1281. *
  1282. * Revision 1.19  2002/04/25 13:14:55  thiessen
  1283. * fix range test
  1284. *
  1285. * Revision 1.18  2002/02/12 13:08:20  thiessen
  1286. * annot description optional
  1287. *
  1288. * Revision 1.17  2002/02/12 12:54:10  thiessen
  1289. * feature legend at bottom; annot only where aligned
  1290. *
  1291. * Revision 1.16  2002/02/08 19:53:17  thiessen
  1292. * add annotation to text/HTML displays
  1293. *
  1294. * Revision 1.15  2001/03/02 01:19:24  thiessen
  1295. * add FASTA output
  1296. *
  1297. * Revision 1.14  2001/02/16 19:18:47  thiessen
  1298. * change color scheme again
  1299. *
  1300. * Revision 1.13  2001/02/15 19:23:44  thiessen
  1301. * add identity coloring
  1302. *
  1303. * Revision 1.12  2001/02/15 18:08:15  thiessen
  1304. * change color scheme
  1305. *
  1306. * Revision 1.11  2001/02/14 16:06:09  thiessen
  1307. * add block and conservation coloring to HTML display
  1308. *
  1309. * Revision 1.10  2001/02/14 03:16:27  thiessen
  1310. * fix seqloc markers and right-justification of left tails
  1311. *
  1312. * Revision 1.9  2001/02/14 00:06:49  thiessen
  1313. * filter out consensus
  1314. *
  1315. * Revision 1.8  2001/01/29 23:55:09  thiessen
  1316. * add AlignmentDisplay verification
  1317. *
  1318. * Revision 1.7  2001/01/29 18:13:33  thiessen
  1319. * split into C-callable library + main
  1320. *
  1321. * Revision 1.6  2001/01/25 00:51:20  thiessen
  1322. * add command-line args; can read asn data from stdin
  1323. *
  1324. * Revision 1.5  2001/01/23 23:17:48  thiessen
  1325. * fix uid link problems
  1326. *
  1327. * Revision 1.4  2001/01/23 20:42:00  thiessen
  1328. * add description
  1329. *
  1330. * Revision 1.3  2001/01/23 17:34:12  thiessen
  1331. * add HTML output
  1332. *
  1333. * Revision 1.2  2001/01/22 15:55:11  thiessen
  1334. * correctly set up ncbi namespacing
  1335. *
  1336. * Revision 1.1  2001/01/22 13:15:23  thiessen
  1337. * initial checkin
  1338. *
  1339. */