/**
 * Author: Mark Larkin
 * 
 * Copyright (c) 2007 Des Higgins, Julie Thompson and Toby Gibson.  
 */
/**
 * Output routines ported from clustalx C code in 'interface.c'.
 *
 * Changes:
 *
 * 18-04-07,Nigel Brown(EMBL): clustalx code used (firsRes, length) convention
 * to define output columns, while the port uses (firsRes, lastRes). Fixed
 * problems in all 7 output routines where the conventions were mixed giving
 * over-long output blocks.
 *
 * 18-7-07: Mark (UCD), made changes to fastaOut, clustalOut, nbrfOut and gdeOut
 *
 * 9-2-08: Paul (UCD), changes to clustalOut to make output the same as 1.83
 *                     added const NAMESWIDTH which is calculated based on the
 *                     length of the longest sequence name and the upper and lower
 *                     limits specified in MAXNAMESTODISPLAY and MINNAMESTODISPLAY
 */

#ifdef HAVE_CONFIG_H
    #include "config.h"
#endif
#include "AlignmentOutput.h"

namespace clustalw
{

AlignmentOutput::AlignmentOutput()
: clusSecStructOffset(9),
  clusSequenceOffset(15+1) // MARK made a change here for test case findRangeValues 10seqs
{
    // NOTE in the old clustal these arrays ended with NULL as a value.
    try
    {
        strongGroup.resize(9);
        strongGroup[0] = string("STA");
        strongGroup[1] = string("NEQK");
        strongGroup[2] = string("NHQK");
        strongGroup[3] = string("NDEQ");
        strongGroup[4] = string("QHRK");
        strongGroup[5] = string("MILV");
        strongGroup[6] = string("MILF");
        strongGroup[7] = string("HY");
        strongGroup[8] = string("FYW");
    
        weakGroup.resize(11);
        weakGroup[0] = string("CSA");
        weakGroup[1] = string("ATV");
        weakGroup[2] = string("SAG");
        weakGroup[3] = string("STNK");
        weakGroup[4] = string("STPA");
        weakGroup[5] = string("SGND");
        weakGroup[6] = string("SNDEQK");
        weakGroup[7] = string("NDEQHK");
        weakGroup[8] = string("NEQHRK");
        weakGroup[9] = string("FVLIM");
        weakGroup[10] = string("HFY");
    }
    catch(const exception& e)
    {
        cerr << "An exception has occured in the contructor of AlignmentOutput.\n"
             << e.what() << "\n";
        exit(1);    
    }
}

/*
 * The function createAlignmentOutput is used to call all the individual functions
 * for the different file types. It is possible to output the alignment in all file
 * types at the same time.
 */
void AlignmentOutput::createAlignmentOutput(Alignment* alignPtr, int firstSeq, int lastSeq)
{
    int length;  
    int firstRes; // starting sequence range - Ramu
    int lastRes; // ending sequence range 
    bool rangeOK;
    
    if((firstSeq <= 0) || (lastSeq < firstSeq))
    {
        utilityObject->error("Cannot produce alignment output."
                " Incorrect call to createAlignmentOutput."
                " firstSeq = %d"
                " lastSeq = %d\n", firstSeq, lastSeq);
        return;
    }
    
    length = 0;
    firstRes = 1;
    lastRes = 0;
    rangeOK = false;
    
    try
    {
        length = alignPtr->getLengthLongestSequence();
        lastRes = length;

        if (userParameters->getRangeFromToSet()) 
        {
            firstRes = userParameters->getRangeFrom();
            lastRes = userParameters->getRangeTo();
            // Check if the numbers are ok.
            if ((firstRes > lastRes) || (firstRes == -1) || (lastRes == -1)) 
            {
                cerr << "seqrange numbers are not set properly, using default....\n";
                firstRes = 1;
                lastRes = length;
            }
            else
            {
                rangeOK = true;
            }
        }
        if (rangeOK && (lastRes > length)) 
        {
            lastRes = length;
            cout << "Seqrange " << lastRes << " is more than the " << length 
                 << "  setting it to " << length << " \n";
        }

        if (userParameters->getMenuFlag())
        {    
            cout << "Consensus length = " << lastRes << " \n";
        }

        outputRegion partToOutput;
        partToOutput._firstSeq = firstSeq;
        partToOutput._lastSeq = lastSeq;
        partToOutput._firstRes = firstRes;
        partToOutput._lastRes = lastRes;

        if(userParameters->getOutputClustal()) 
        {
            if(clustalOutFile.get() != 0 && clustalOutFile->is_open())
            {
                clustalOut(alignPtr, partToOutput);
                if(clustalOutFile->is_open())
                {
                    clustalOutFile->close();
                }
                utilityObject->info("CLUSTAL-Alignment file created  [%s]\n",
                                    clustalOutName.c_str());
            }
        }
        if(userParameters->getOutputNbrf())  
        {
            if(nbrfOutFile.get() != 0 && nbrfOutFile->is_open())
            {            
                nbrfOut(alignPtr, partToOutput);
                if(nbrfOutFile->is_open())
                {
                    nbrfOutFile->close();
                }
                utilityObject->info("NBRF/PIR-Alignment file created  [%s]\n",
                                    nbrfOutName.c_str());                
            }
        }
        if(userParameters->getOutputGCG())  
        {
            if(gcgOutFile.get() != 0 && gcgOutFile->is_open())
            {
                gcgOut(alignPtr, partToOutput);
                if(gcgOutFile->is_open())
                {
                    gcgOutFile->close();
                }
                utilityObject->info("GCG-Alignment file created      [%s]\n",
                                    gcgOutName.c_str());                
            }
        }
        if(userParameters->getOutputPhylip())  
        {
            if(phylipOutFile.get() != 0 && phylipOutFile->is_open())
            {
                phylipOut(alignPtr, partToOutput);
                if(phylipOutFile->is_open())
                {
                    phylipOutFile->close();
                }
                utilityObject->info("PHYLIP-Alignment file created   [%s]\n",
                                    phylipOutName.c_str());                
            }
        }
        
        if(userParameters->getOutputGde())  
        {
            if(gdeOutFile.get() != 0 && gdeOutFile->is_open())
            {
                gdeOut(alignPtr, partToOutput);
                if(gdeOutFile->is_open())
                {
                    gdeOutFile->close();
                }
                utilityObject->info("GDE-Alignment file created      [%s]\n",
                                    gdeOutName.c_str());                
            }
        }
        
        if(userParameters->getOutputNexus())  
        {
            if(nexusOutFile.get() != 0 && nexusOutFile->is_open())
            {
                nexusOut(alignPtr, partToOutput);
                if(nexusOutFile->is_open())
                {
                    nexusOutFile->close();
                }
                utilityObject->info("NEXUS-Alignment file created    [%s]\n",
                                    nexusOutName.c_str());                
            }
        }

        if(userParameters->getOutputFasta())  
        {
            if(fastaOutFile.get() != 0 && fastaOutFile->is_open())
            {
                fastaOut(alignPtr, partToOutput);
                if(fastaOutFile->is_open())
                {
                    fastaOutFile->close();
                }
                utilityObject->info("Fasta-Alignment file created    [%s]\n",
                                    fastaOutName.c_str());                
            }
        }
        
        // Output the alignment to the screen if it is required.
        if (userParameters->getShowAlign() && userParameters->getMenuFlag()) 
        {
            showAlign();
        }
    }
    catch(const exception& e)
    {
        cerr << "An exception has occured in the createAlignmentOutput function.\n"
             << e.what() << "\n";
        exit(1);    
    }
}

bool AlignmentOutput::QTOpenFilesForOutput(AlignmentFileNames fileNames)
{
    if(!userParameters->getOutputClustal() && 
       !userParameters->getOutputNbrf() && !userParameters->getOutputGCG() &&
       !userParameters->getOutputPhylip() && !userParameters->getOutputGde() &&
       !userParameters->getOutputNexus() && !userParameters->getOutputFasta()) 
    {
         utilityObject->error("You must select an alignment output format\n");
         return false;
    }
        
    if(fileNames.clustalFile == "" && fileNames.fastaFile == "" && fileNames.gcgFile == ""
       && fileNames.gdeFile == "" && fileNames.nexusFile == "" && fileNames.nrbfFile == ""
       && fileNames.phylipFile == "")
    {
        utilityObject->error("No names for output files. Cannot output alignment.\n");
        return false;
    }
    
    if(fileNames.clustalFile != "")
    {
        clustalOutName = fileNames.clustalFile;
        if(!openExplicitFile(clustalOutFile, clustalOutName))
        { 
            return false;
        }
    }  
    if(fileNames.fastaFile != "")
    {
        fastaOutName = fileNames.fastaFile;
        if(!openExplicitFile(fastaOutFile, fastaOutName))
        { 
            return false;
        }
    }
    if(fileNames.gcgFile != "")
    {
        gcgOutName = fileNames.gcgFile;
        if(!openExplicitFile(gcgOutFile, gcgOutName))
        { 
            return false;
        }
    }
    if(fileNames.gdeFile != "")
    {
        gdeOutName = fileNames.gdeFile;
        if(!openExplicitFile(gdeOutFile, gdeOutName))
        { 
            return false;
        }
    }
    if(fileNames.nexusFile != "")
    {
        nexusOutName = fileNames.nexusFile;
        if(!openExplicitFile(nexusOutFile, nexusOutName))
        { 
            return false;
        }
    }
    if(fileNames.nrbfFile != "")
    {
        nbrfOutName = fileNames.nrbfFile;
        if(!openExplicitFile(nbrfOutFile, nbrfOutName))
        { 
            return false;
        }
    }
    if(fileNames.phylipFile != "")
    {
        phylipOutName = fileNames.phylipFile;
        if(!openExplicitFile(phylipOutFile, phylipOutName))
        { 
            return false;
        }
    } 
    return true;                         
}
/*
 * The function openAlignmentOutput opens a file for output. It returns true if it
 * has been successful and false if it has not been successful
 */
bool AlignmentOutput::openAlignmentOutput(string path)
{
    if(!userParameters->getOutputClustal() && 
       !userParameters->getOutputNbrf() && !userParameters->getOutputGCG() &&
       !userParameters->getOutputPhylip() && !userParameters->getOutputGde() &&
       !userParameters->getOutputNexus() && !userParameters->getOutputFasta()) 
    {
         utilityObject->error("You must select an alignment output format\n");
         return false;
    }
    string _fileNameToOutput = path;
    if(_fileNameToOutput == "")
    {
      /*_fileNameToOutput = userParameters->getSeqName();*/
      /* BUG 166 , file extension, FS, 2009-01-26 */
      utilityObject->getPath(userParameters->getSeqName(), &_fileNameToOutput);
    }
    
    if(userParameters->getOutputClustal()) 
    {
        if (userParameters->getOutfileName() != "") 
        {
            clustalOutName = userParameters->getOutfileName();
            if(!openExplicitFile(clustalOutFile, clustalOutName))
            { 
                return false;
            }
        }
        else 
        {
            clustalOutName = openOutputFile(clustalOutFile, 
                             "\nEnter a name for the CLUSTAL output file ",
                             _fileNameToOutput, "aln");
                             
            if(clustalOutName == "")
            {
                return false; // We have not been successful.
            }
        }
    }

    if(userParameters->getOutputNbrf()) 
    {
        if (userParameters->getOutfileName() != "") 
        {
            nbrfOutName = userParameters->getOutfileName();
            if(!openExplicitFile(nbrfOutFile, nbrfOutName))
            { 
                return false;
            }            
        }
        else
        {
            nbrfOutName = openOutputFile(nbrfOutFile,
                      "\nEnter a name for the NBRF/PIR output file", _fileNameToOutput,
                       "pir");
                      
            if(nbrfOutName == "")
            {
                return false; // We have not been successful.
            }
        }
    }

    if(userParameters->getOutputGCG()) 
    {
        if (userParameters->getOutfileName() != "") 
        {
            gcgOutName = userParameters->getOutfileName();
            if(!openExplicitFile(gcgOutFile, gcgOutName))
            { 
                return false;
            }   
        }
        else
        {
            gcgOutName = openOutputFile(gcgOutFile,
                      "\nEnter a name for the GCG output file", _fileNameToOutput, "msf");
                      
            if(gcgOutName == "")
            {
                return false; // We have not been successful.
            }
        }
    }

    if(userParameters->getOutputPhylip()) 
    {
        if (userParameters->getOutfileName() != "") 
        {
            phylipOutName = userParameters->getOutfileName();
            if(!openExplicitFile(phylipOutFile, phylipOutName))
            { 
                return false;
            }
        }
        else
        {
            phylipOutName = openOutputFile(phylipOutFile,
                      "\nEnter a name for the PHYLIP output file", _fileNameToOutput, "phy");
                      
            if(phylipOutName == "")
            {
                return false; // We have not been successful.
            }
        }
    }

    if(userParameters->getOutputGde()) 
    {
        if (userParameters->getOutfileName() != "") 
        {
            gdeOutName = userParameters->getOutfileName();
            if(!openExplicitFile(gdeOutFile, gdeOutName))
            { 
                return false;
            }
        }
        else
        {
            gdeOutName = openOutputFile(gdeOutFile, 
                             "\nEnter a name for the GDE output file     ",
                             _fileNameToOutput, "gde");
                 
            if(gdeOutName == "")
            {
                return false; // We have not been successful.
            }                       
        }
    }

    if(userParameters->getOutputNexus()) 
    {
        if (userParameters->getOutfileName() != "") 
        {
            nexusOutName = userParameters->getOutfileName();
            if(!openExplicitFile(nexusOutFile, nexusOutName))
            { 
                return false;
            }
        }
        else
        {
            nexusOutName = openOutputFile(nexusOutFile, 
                             "\nEnter a name for the NEXUS output file   ",
                             _fileNameToOutput, "nxs");
                 
            if(nexusOutName == "")
            {
                return false; // We have not been successful.
            }                       
        }
    }
   
    if(userParameters->getOutputFasta()) 
    {
        if (userParameters->getOutfileName() != "") 
        {
            fastaOutName = userParameters->getOutfileName();
            if(!openExplicitFile(fastaOutFile, fastaOutName))
            { 
                return false;
            }
        }
        else
        {
            fastaOutName = openOutputFile(fastaOutFile, 
                             "\nEnter a name for the Fasta output file ",
                             _fileNameToOutput, "fasta");
                 
            if(fastaOutName == "")
            {
                return false; // We have not been successful.
            }
        }
    }
  
    return true;

}


/*
 * The function openExplicitFile is used to open a file when we have the name already.
 * It returns true if the file has been opened correctly, false otherwise.
 */
bool AlignmentOutput::openExplicitFile(auto_ptr<ofstream>& outFile, string fileName)
{
    if (fileName == "") 
    {
        cerr << "Bad output file [" << fileName << "]\n";
        utilityObject->error("Bad output file [%s]\n", fileName.c_str());
        return false;
    }
    
    outFile.reset(new ofstream(fileName.c_str(), ofstream::trunc));
    if(!outFile->is_open()) 
    {
        utilityObject->error("Cannot open output file [%s]\n", fileName.c_str());
        return false;
    }
    return true;
}

/*
 * The function openOutputFile is used to open a file when we dont have the name yet.
 * This function returns the name if it has been successful. If it has not been 
 * successful it returns a blank string ""
 */
string AlignmentOutput::openOutputFile(auto_ptr<ofstream>& outFile, string prompt, 
                                       string path, string fileExtension)
{
    string temp;
    string _fileName; // Will return this name.
    string message;
    _fileName = path + fileExtension;

    if(_fileName.compare(userParameters->getSeqName()) == 0) 
    {
        cout << "Output file name is the same as input file.\n";
        if (userParameters->getMenuFlag()) 
        {
            message = "\n\nEnter new name to avoid overwriting  [" + _fileName + "]";

            utilityObject->getStr(message, temp);
            if(temp != "")
            {
                _fileName = temp;
            }
        }
    }
    else if (userParameters->getMenuFlag()) 
    {

        message = prompt + " [" + _fileName + "]";
        utilityObject->getStr(message, temp);
        if(temp != "")
        {
            _fileName = temp;
        }
    }

    outFile.reset(new ofstream(_fileName.c_str(), ofstream::trunc));
    if(!outFile->is_open()) 
    {
        utilityObject->error("Cannot open output file [%s]\n", _fileName.c_str());
        return "";
    }    
    return _fileName;
}

/*
 * fastaOut: output the alignment in FASTA format
 */
void AlignmentOutput::fastaOut(Alignment* alignPtr, outputRegion partToOutput)
{
    char residue;
    int val;
    int i, ii;
    int j, slen;    
    int lineLength;
    int firstRes = partToOutput._firstRes;
    int lastRes = partToOutput._lastRes;
    int firstSeq = partToOutput._firstSeq;
    int lastSeq = partToOutput._lastSeq;
    rangeNum rnum;  
    cout << "firstres = " << firstRes << " lastres = " << lastRes << "\n";
    try
    {
        const SeqArray* alignment = alignPtr->getSeqArray(); //NOTE june29
        vector<char> sequence;
        sequence.assign(lastRes + 1, '0');     
    
        lineLength = PAGEWIDTH - alignPtr->getMaxNames();
                
        lineLength = lineLength - lineLength % 10; // round to a multiple of 10
        if (lineLength > LINELENGTH || lineLength <= 0) // Mark 18-7-07
        {    
            lineLength = LINELENGTH;
        }


        for(ii = firstSeq; ii <= lastSeq; ii++) 
        {
            i = alignPtr->getOutputIndex(ii - 1);
            slen = 0;
            //for(j = firstRes; j < firstRes + lastRes; j++) //- nige
            for(j = firstRes; j <= lastRes; j++) //- nige
            {
                if (j <= alignPtr->getSeqLength(i))
                {
                    //val = alignPtr.getResidue(i, j);
                    val = (*alignment)[i][j]; // NOTE june29
                }
                else
                { 
                    val = -3;
                }
                if((val == -3) || (val == 253))
                {
                    break;
                }
                else if((val < 0) || (val > userParameters->getMaxAA())) 
                {
                    residue = '-';
                }
                else 
                {
                    residue = userParameters->getAminoAcidCode(val);
                }
                if (userParameters->getLowercase())
                {    
                    sequence[j - firstRes] = residue;
                }
                else
                {
                    sequence[j - firstRes] = residue;
                }
                slen++;
            }
            (*fastaOutFile) << ">" << nameonly(alignPtr->getName(i));
            if(userParameters->getSeqRange()) 
            {
                findRangeValues(alignPtr, &rnum, firstRes, lastRes, ii);
                (*fastaOutFile) << "/" << rnum.start << "-" << rnum.end;
                
            }
            (*fastaOutFile) << "\n";
            for(j = 1; j <= slen; j++) 
            {
                (*fastaOutFile) << sequence[j-1];
                if((j % lineLength == 0) || (j == slen))
                { 
                    (*fastaOutFile) << "\n";
                }
            }
        }
        fastaOutFile->close();
        cout << "FASTA file created!\n";
    }
    catch(const bad_alloc& e)
    {
        fastaOutFile->close();
        cerr << "A bad_alloc exception has occured in the fastaOut function.\n"
             << e.what() << "\n";
        exit(1);
    }
    catch(VectorOutOfRange e)
    {
        fastaOutFile->close();
        cerr << "An exception has occured in the fastaOut function.\n"
             << e.what() << "\n";
        exit(1);    
    }
    catch(...)
    {
        fastaOutFile->close();
        cerr << "An exception has occured in the fastaOut function.\n";
        exit(1);    
    }

}

/* 
 * gcgOut: output the alignment in gcg file format
 */
void AlignmentOutput::gcgOut(Alignment* alignPtr, outputRegion partToOutput)
{
    char residue;
    int val;
    int i, ii, chunks, block;
    int j, k, pos1, pos2;    
    long grandChecksum;
 
    try
    {  
        int firstRes = partToOutput._firstRes;
        int lastRes = partToOutput._lastRes;
        int firstSeq = partToOutput._firstSeq;
        int lastSeq = partToOutput._lastSeq;    
        rangeNum rnum;
        vector<char> sequence;
        vector<int> allChecks;
        int _maxLength = alignPtr->getMaxAlnLength();
        sequence.assign(_maxLength + 1, '0');
        allChecks.assign(lastSeq + 1, 0);
        const SeqArray* alignment = alignPtr->getSeqArray();
        
        for(i = firstSeq; i <= lastSeq; i++) 
        {
            //for(j = firstRes; j <= firstRes + lastRes - 1; j++) //- nige
            for(j = firstRes; j <= lastRes; j++) //- nige
            {
                if (j <= alignPtr->getSeqLength(i))
                {
                    //val = alignPtr.getResidue(i, j);
                    val = (*alignment)[i][j]; // NOTE june29
                }
                else
                { 
                    val = -3;
                }
                if((val == -3) || (val == 253)) 
                {
                    break;
                }
                else if((val < 0) || (val > userParameters->getMaxAA()))
                {
                    residue = '.';
                }
                else 
                {
                    residue = userParameters->getAminoAcidCode(val);
                }
                sequence[j - firstRes + 1] = residue;
            }
            // pad any short sequences with gaps, to make all sequences the same length 
            for(; j <= firstRes + lastRes - 1; j++)
            { 
                sequence[j - firstRes + 1] = '.';
            }
            allChecks[i] = SeqGCGCheckSum(&sequence, lastRes);
        }    
        int _index;
        grandChecksum = 0;
        for(i = 1; i <= alignPtr->getNumSeqs(); i++)
        {
            _index = alignPtr->getOutputIndex(i - 1);
            grandChecksum += allChecks[_index];
        }
        grandChecksum = grandChecksum % 10000;
        
        (*gcgOutFile) << "PileUp\n\n";
        (*gcgOutFile) << "\n\n   MSF:" << setw(5) << lastRes << "  Type: ";
                        
        if(userParameters->getDNAFlag())
        {
            (*gcgOutFile) << "N";
        }
        else
        {
            (*gcgOutFile) << "P";
        }

        (*gcgOutFile) << "    Check:" << setw(6) << grandChecksum << "   .. \n\n";
        float _seqWeight;

        int length = lastRes - firstRes + 1; //- nige

        for(ii = firstSeq; ii <= lastSeq; ii++)  
        {
            i = alignPtr->getOutputIndex(ii - 1);
            _seqWeight = (alignPtr->getSeqWeight(i - 1) * 100) / (float) INT_SCALE_FACTOR;

            (*gcgOutFile) << " Name: " << alignPtr->getName(i) << " oo  Len:"
                        //<< setw(5) << lastRes << "  Check:" << setw(6) << allChecks[i]  //- nige
                          << setw(5) << length << "  Check:" << setw(6) << allChecks[i] 
                          << "  Weight:  " << fixed << setprecision(1) << _seqWeight << "\n";
        }  
        (*gcgOutFile) << "\n//\n";
        
        chunks = length / GCG_LINELENGTH; //- nige
        if(length % GCG_LINELENGTH != 0)  //- nige
        {
            ++chunks;
        }
  
        for(block = 1; block <= chunks; block++) 
        {
            (*gcgOutFile) << "\n\n";
            pos1 = ((block - 1) * GCG_LINELENGTH) + 1;
            pos2 = (length < pos1 + GCG_LINELENGTH - 1)? length : pos1 + GCG_LINELENGTH - 1;//- nige
            
            for(ii = firstSeq; ii <= lastSeq; ii++) 
            {
                i = alignPtr->getOutputIndex(ii - 1);
                if (!userParameters->getSeqRange()) 
                {
                    (*gcgOutFile) << "\n" << setw(alignPtr->getMaxNames() + 5) << left
                                  << alignPtr->getName(i) << " ";
                }
                else 
                {
                    findRangeValues(alignPtr, &rnum, firstRes, lastRes, ii);
                    std::stringstream ss;
                    std::stringstream ss2;
                    string rangeStart;
                    string rangeEnd;
                    ss << rnum.start;
                    ss >> rangeStart;
                    ss2 << rnum.end;
                    ss2 >> rangeEnd;
                    string nameAndRange = nameonly(alignPtr->getName(i)) + "/" + rangeStart;
                    nameAndRange += "-" + rangeEnd;
                    (*gcgOutFile) << "\n" << setw(alignPtr->getMaxNames() + 15) << left
                                  << nameAndRange;
                }
                for(j = pos1, k = 1; j <= pos2; j++, k++) 
                {
            
                    //JULIE -
                    //check for int sequences - pad out with '.' characters to end 
            
                    if (j + firstRes - 1 <= alignPtr->getSeqLength(i))
                    {
                        //val = alignPtr.getResidue(i, j + firstRes - 1);
                        val = (*alignment)[i][j + firstRes - 1]; // NOTE june29
                    }
                    else
                    {
                        val = -3;
                    }
                    if((val == -3) || (val == 253))
                    {
                        residue = '.';
                    }
                    else if((val < 0) || (val > userParameters->getMaxAA()))
                    {
                        residue = '.';
                    }
                    else 
                    {
                        residue = userParameters->getAminoAcidCode(val);
                    }
                    (*gcgOutFile) << residue;
                    if(j % 10 == 0)
                    { 
                        (*gcgOutFile) << " ";
                    }
                }
            }
        }
        (*gcgOutFile) << "\n\n";
        gcgOutFile->close();   
    }
    catch(const bad_alloc& e)
    {
        gcgOutFile->close();
        cerr << "A bad_alloc exception has occured in the gcgOut function.\n"
             << e.what() << "\n";
        exit(1);
    }
    catch(VectorOutOfRange e)
    {
        gcgOutFile->close();
        cerr << "An exception has occured in the gcgOut function.\n"
             << e.what() << "\n";
        exit(1);    
    }
    catch(...)
    {
        gcgOutFile->close();
        cerr << "An exception has occured in the gcgOut function.\n";
        exit(1);    
    }    
}

void AlignmentOutput::phylipOut(Alignment* alignPtr, outputRegion partToOutput)
{
    int firstRes = partToOutput._firstRes;
    int lastRes = partToOutput._lastRes;
    int firstSeq = partToOutput._firstSeq;
    int lastSeq = partToOutput._lastSeq;    
    try
    {
        char residue;
        int val;
        int i, ii, chunks, block;    
        int j, k, pos1, pos2;    
        int nameLen;
        bool warn;
        vector<string> _seqNames;
        const SeqArray* alignment = alignPtr->getSeqArray();
        rangeNum rnum;

        // Push on a blank string. This is because the index starts at 1 for the seqs etc..
        _seqNames.push_back("");

        nameLen = 0;
    
        for(i = firstSeq; i <= lastSeq; i++)  
        {
            _seqNames.push_back(alignPtr->getName(i));
            ii = _seqNames.size();
            if(nameLen < ii) 
            {
                nameLen = ii;
            }
        }
        if(nameLen > 10) 
        {
            warn = false;
            for(i = 0; i < (int)_seqNames.size() - 1; i++)  
            {
                for(j = i + 1; j < (int)_seqNames.size(); j++) 
                {
                    if (_seqNames[i].compare(_seqNames[j]) == 0) 
                    {
                        warn = true;
                    }
                }
            }
            if(warn)
            {
                utilityObject->warning("Truncating sequence names to 10 characters for PHYLIP output.\nNames in the PHYLIP format file are NOT unambiguous.");
            }
            else
            {
                utilityObject->warning("Truncating sequence names to 10 characters for PHYLIP output.");
            }
        }
  
        int length = lastRes - firstRes + 1; //- nige

        chunks = length / GCG_LINELENGTH; //- nige
        if(length % GCG_LINELENGTH != 0)  //- nige
        {
            ++chunks;
        }
  
        (*phylipOutFile) << setw(6) << alignPtr->getNumSeqs() << " "
                         << setw(6) << length; //-nige
        
        for(block = 1; block <= chunks; block++) 
        {
            pos1 = ((block - 1) * GCG_LINELENGTH) + 1;
            //pos2 = (lastRes < pos1 + GCG_LINELENGTH - 1)? lastRes : pos1 + GCG_LINELENGTH - 1; //- nige
            pos2 = (length < pos1 + GCG_LINELENGTH - 1)? length : pos1 + GCG_LINELENGTH - 1;
            for(ii = firstSeq; ii <= lastSeq; ii++)  
            {
                i = alignPtr->getOutputIndex(ii - 1);
                if(block == 1)  
                {
                    if(!userParameters->getSeqRange()) 
                    {
                        string name = alignPtr->getName(i);
                        (*phylipOutFile) << "\n" << setw(10) << left 
                                             << name.substr(0,10) << " ";
                    }
                    else
                    {
                        findRangeValues(alignPtr, &rnum, firstRes, lastRes, ii);
                        std::stringstream ss;
                        std::stringstream ss2;
                        string rangeStart;
                        string rangeEnd;
                        ss << rnum.start;
                        ss >> rangeStart;
                        ss2 << rnum.end;
                        ss2 >> rangeEnd;
                        string nameAndRange = nameonly(alignPtr->getName(i)) + "/";
                        nameAndRange += rangeStart + "-" + rangeEnd;
                        (*phylipOutFile) << "\n" << setw(alignPtr->getMaxNames() + 15)
                                         << left
                                         << nameAndRange;
                    }
                }
                else
                {
                    (*phylipOutFile) << "\n           ";
                }
                for(j = pos1, k = 1; j <= pos2; j++, k++) 
                {
                    if (j + firstRes - 1 <= alignPtr->getSeqLength(i))
                    {
                        //val = alignPtr.getResidue(i, j + firstRes - 1);
                        val = (*alignment)[i][j + firstRes - 1]; // NOTE june29
                    }
                    else
                    {
                        val = -3;
                    }
                    if((val == -3) || (val == 253))
                    {
                        break;
                    }
                    else if((val < 0) || (val > userParameters->getMaxAA()))
                    {
                        residue = '-';
                    }
                    else 
                    {
                        residue = userParameters->getAminoAcidCode(val);
                    }
                    (*phylipOutFile) << residue;
                    if(j % 10 == 0) 
                    {
                        (*phylipOutFile) << " ";
                    }
                }
            }
            (*phylipOutFile) << "\n";
        }  
    }
    catch(const bad_alloc& e)
    {
        phylipOutFile->close();
        cerr << "A bad_alloc exception has occured in the phylipOut function.\n"
             << e.what() << "\n";
        exit(1);
    }
    catch(VectorOutOfRange e)
    {
        phylipOutFile->close();
        cerr << "An exception has occured in the phylipOut function.\n"
             << e.what() << "\n";
        exit(1);    
    }
    catch(...)
    {
        phylipOutFile->close();
        cerr << "An exception has occured in the phylipOut function.\n";
        exit(1);    
    }      
}

void AlignmentOutput::nexusOut(Alignment* alignPtr, outputRegion partToOutput)
{
    int firstRes = partToOutput._firstRes;
    int lastRes = partToOutput._lastRes;
    int firstSeq = partToOutput._firstSeq;
    int lastSeq = partToOutput._lastSeq;    
    try
    {
        char residue;
        int val;
        int i, ii, chunks, block;    
        int j, k, pos1, pos2;    
        const SeqArray* alignment = alignPtr->getSeqArray();
        rangeNum rnum;

        int length = lastRes - firstRes + 1; //- nige

        chunks = length / GCG_LINELENGTH; //- nige
        if(length % GCG_LINELENGTH != 0) //- nige
        {
            ++chunks;
        }
  
        (*nexusOutFile) << "#NEXUS\n";
        (*nexusOutFile) << "BEGIN DATA;\n";
        (*nexusOutFile) << "dimensions ntax=" << alignPtr->getNumSeqs() 
                        << " nchar=" << length << ";\n"; //- nige
        (*nexusOutFile) << "format missing=?\n";
        (*nexusOutFile) << "symbols=\"";
        
        for(i = 0; i <= userParameters->getMaxAA(); i++)
        {
            (*nexusOutFile) << userParameters->getAminoAcidCode(i);
        }
        (*nexusOutFile) << "\"\n";
        (*nexusOutFile) << "interleave datatype=";
        bool _dnaFlag = userParameters->getDNAFlag();
        string _type = _dnaFlag ? "DNA " : "PROTEIN ";
        (*nexusOutFile) << _type;
        (*nexusOutFile) << "gap= -;\n";
        (*nexusOutFile) << "\nmatrix";
  
        for(block = 1; block <= chunks; block++) 
        {
            pos1 = ((block - 1) * GCG_LINELENGTH) + 1;
            pos2 = (length < pos1 + GCG_LINELENGTH - 1)? length : pos1 + GCG_LINELENGTH - 1; //- nige
            for(ii = firstSeq; ii <= lastSeq; ii++)  
            {
                i = alignPtr->getOutputIndex(ii - 1);
                if (!userParameters->getSeqRange()) 
                {
                    (*nexusOutFile) << "\n" << setw(alignPtr->getMaxNames() + 1) << left
                                    << alignPtr->getName(i) << " ";
                }
                else 
                {
                    findRangeValues(alignPtr, &rnum, firstRes, lastRes, ii);
                    
                    std::stringstream ss;
                    std::stringstream ss2;
                    string rangeStart;
                    string rangeEnd;
                    ss << rnum.start;
                    ss >> rangeStart;
                    ss2 << rnum.end;
                    ss2 >> rangeEnd;
                    string nameAndRange = nameonly(alignPtr->getName(i)) + "/";
                    nameAndRange += rangeStart + "-" + rangeEnd;
                    (*nexusOutFile) << "\n" << setw(alignPtr->getMaxNames() + 15)
                                    << left
                                    << nameAndRange;                    
                }
                for(j = pos1, k = 1; j <= pos2; j++, k++) 
                {
                    if (j + firstRes - 1 <= alignPtr->getSeqLength(i))
                    {
                        //val = alignPtr.getResidue(i, j + firstRes - 1);
                        val = (*alignment)[i][j + firstRes - 1]; // NOTE june29
                    }
                    else
                    {
                        val = -3;
                    }
                    if((val == -3) || (val == 253))
                    {
                        break;
                    }
                    else if((val < 0) || (val > userParameters->getMaxAA()))
                    {
                        residue = '-';
                    }
                    else 
                    {
                        residue = userParameters->getAminoAcidCode(val);
                    }
                    (*nexusOutFile) << residue;
                }
            }
            (*nexusOutFile) << "\n";
        }
        (*nexusOutFile) << ";\nend;\n";   
    }
    catch(const bad_alloc& e)
    {
        nexusOutFile->close();
        cerr << "A bad_alloc exception has occured in the nexusOut function.\n"
             << e.what() << "\n";
        exit(1);
    }
    catch(VectorOutOfRange e)
    {
        nexusOutFile->close();
        cerr << "An exception has occured in the nexusOut function.\n"
             << e.what() << "\n";
        exit(1);    
    }
    catch(...)
    {
        nexusOutFile->close();
        cerr << "An exception has occured in the nexusOut function.\n";
        exit(1);    
    }      
}

/*
 * gdeOut: print out the alignment in gde file format.
 */
void AlignmentOutput::gdeOut(Alignment* alignPtr, outputRegion partToOutput)
{
    int firstRes = partToOutput._firstRes;
    int lastRes = partToOutput._lastRes;
    int firstSeq = partToOutput._firstSeq;
    int lastSeq = partToOutput._lastSeq;    

    int length = lastRes - firstRes + 1; //- nige

    try
    {
        char residue;
        int val;
        vector<char> _ssMask1, _ssMask2, seq;
        int i, ii;
        int j, slen;    
        int lineLength;
        const SeqArray* alignment = alignPtr->getSeqArray();
        rangeNum rnum;

        seq.assign(alignPtr->getMaxAlnLength() + 1, '0');
        
        // decide the line length for this alignment - maximum is LINELENGTH 
        lineLength = PAGEWIDTH - alignPtr->getMaxNames();
        lineLength = lineLength - lineLength % 10; // round to a multiple of 10
        if (lineLength > LINELENGTH || lineLength <= 0) 
        {
            lineLength = LINELENGTH;
        }
  
        // If we are using the secondary structure info from profile 1, set it up
        if (userParameters->getStructPenalties1() == SECST && 
            userParameters->getUseSS1() == true) 
        {
            int _lengthSeq1 = alignPtr->getSeqLength(1);
            vector<char>* _secStructMask1 = alignPtr->getSecStructMask1();
            _ssMask1.assign(_lengthSeq1 + 10, 0);
            
            for (i = 0;i < _lengthSeq1;i++)
            {
                _ssMask1[i] = _secStructMask1->at(i);
            }
            printSecStructMask(_lengthSeq1, _secStructMask1, &_ssMask1);
        }
        
        // If we are using the secondary structure info from profile 2, set it up
        if (userParameters->getStructPenalties2() == SECST && 
            userParameters->getUseSS2() == true) 
        {
            int indexProfile2FirstSeq = alignPtr->getProfile1NumSeqs() + 1;
            int lengthSeqProfile2 = alignPtr->getSeqLength(indexProfile2FirstSeq);
            _ssMask2.assign(lengthSeqProfile2 + 10, 0);
            vector<char>* _secStructMask2 = alignPtr->getSecStructMask2();
            
            for (i=0;i < lengthSeqProfile2;i++)
            {
                _ssMask2[i] = _secStructMask2->at(i);
            }
            printSecStructMask(lengthSeqProfile2, _secStructMask2, &_ssMask2);  
        }

        bool _dnaFlag = userParameters->getDNAFlag();
        string _prefix = _dnaFlag ? "#" : "%";
        
        for(ii = firstSeq; ii <= lastSeq; ii++) 
        {
            i = alignPtr->getOutputIndex(ii - 1);
            
            (*gdeOutFile) << _prefix;
            if(!userParameters->getSeqRange()) 
            {
                (*gdeOutFile) << alignPtr->getName(i) << "\n";
            }
            else 
            {
                findRangeValues(alignPtr, &rnum, firstRes, lastRes, ii);
                (*gdeOutFile) << nameonly(alignPtr->getName(i)) << "/"
                              << rnum.start << "-" << rnum.end << "\n";
            }
            slen = 0;

            //for(j = firstRes; j < firstRes + lastRes; j++) //- nige
            for(j = firstRes; j <= lastRes; j++) //- nige
            {
                if (j <= alignPtr->getSeqLength(i))
                {
                    //val = alignPtr.getResidue(i, j);
                    val = (*alignment)[i][j]; // NOTE june29
                }
                else
                { 
                    val = -3;
                }
                if((val == -3) || (val == 253))
                {
                    break;
                }
                else if((val < 0) || (val > userParameters->getMaxAA()))
                {
                    residue = '-';
                }
                else 
                {
                    residue = userParameters->getAminoAcidCode(val);
                }
                if (userParameters->getLowercase())
                {
                    seq[j - firstRes] = (char)tolower((int)residue);
                }
                else
                {
                    seq[j - firstRes] = residue;
                }
                slen++;
            }
            for(j = 1; j <= slen; j++) 
            {
                (*gdeOutFile) << seq[j-1];
                if((j % lineLength == 0) || (j == slen)) 
                {
                    (*gdeOutFile) << "\n";
                }
            }
        }
  
        if (userParameters->getOutputStructPenalties() == 0 || 
            userParameters->getOutputStructPenalties() == 2) 
        {
            if (userParameters->getStructPenalties1() == SECST && 
                userParameters->getUseSS1() == true) 
            {
                (*gdeOutFile) << "\"SS_" << setw(alignPtr->getMaxNames()) << left
                              << alignPtr->getSecStructName1() << "\n";
                              
                //for(i = firstRes; i < firstRes + lastRes; i++) //- nige
                for(i = firstRes; i <= lastRes; i++) //- nige
                {
                    val = _ssMask1[i - 1];
                    if (val == userParameters->getGapPos1() || 
                        val == userParameters->getGapPos2())
                    {
                        seq[i - firstRes] = '-';
                    }
                    else
                    {
                        seq[i - firstRes] = val;
                    }
                }

                for(i = 1; i <= lastRes; i++) 
                {
                    (*gdeOutFile) << seq[i-1];
                    if((i % lineLength == 0) || (i == lastRes)) 
                    {
                        (*gdeOutFile) << "\n";
                    }
                }
            }
    
            if (userParameters->getStructPenalties2() == SECST && 
                userParameters->getUseSS2() == true) 
            {
                (*gdeOutFile) << "\"SS_" << setw(alignPtr->getMaxNames()) << left
                              << alignPtr->getSecStructName2() << "\n";
                                              
                //for(i = firstRes; i < firstRes + lastRes; i++) //- nige
                for(i = firstRes; i <= lastRes; i++) //- nige
                {
                    val = _ssMask2[i - 1];
                    if (val == userParameters->getGapPos1() || 
                        val == userParameters->getGapPos2())
                    {
                        seq[i - firstRes] = '-';
                    }
                    else
                    {
                        seq[i - firstRes] = val;
                    }
                }

                //for(i = 1; i <= lastRes; i++) //- nige
                for(i = 1; i <= length; i++) //- nige
                {
                    (*gdeOutFile) << seq[i - 1];
                    //if((i % lineLength == 0) || (i == lastRes)) //- nige
                    if((i % lineLength == 0) || (i == length)) //- nige
                    {
                        (*gdeOutFile) << "\n";
                    }
                }
            }
        }
        if (userParameters->getOutputStructPenalties() == 1 || 
            userParameters->getOutputStructPenalties() == 2) 
        {
            if (userParameters->getStructPenalties1() != NONE && 
                userParameters->getUseSS1() == true) 
            {
                (*gdeOutFile) << "\"GM_" << setw(alignPtr->getMaxNames()) << left
                              << alignPtr->getSecStructName1() << "\n";
                
                vector<char>* _gapPenaltyMask1 = alignPtr->getGapPenaltyMask1();          
                //for(i = firstRes; i < firstRes + lastRes; i++) //- nige
                for(i = firstRes; i <= lastRes; i++) //- nige
                {
                    val = _gapPenaltyMask1->at(i - 1);
                    if (val == userParameters->getGapPos1() || 
                        val == userParameters->getGapPos2())
                    {
                        seq[i - firstRes] = '-';
                    }
                    else
                    {
                        seq[i - firstRes] = val;
                    }
                }

                //for(i = 1; i <= lastRes; i++) //- nige
                for(i = 1; i <= length; i++) //- nige
                {
                    (*gdeOutFile) << seq[i - 1];
                    //if((i % lineLength == 0) || (i == lastRes)) //- nige
                    if((i % lineLength == 0) || (i == length)) //- nige
                    {
                        (*gdeOutFile) << "\n";
                    }
                }
            }
            if (userParameters->getStructPenalties2() != NONE && 
                userParameters->getUseSS2() == true) 
            {
                (*gdeOutFile) << "\"GM_" << setw(alignPtr->getMaxNames()) << left
                              << alignPtr->getSecStructName2() << "\n";
                
                vector<char>* _gapPenaltyMask2 = alignPtr->getGapPenaltyMask2();    
                //for(i = firstRes; i < firstRes + lastRes; i++) //- nige
                for(i = firstRes; i < length; i++) //- nige
                {
                    val = _gapPenaltyMask2->at(i-1);
                    if (val == userParameters->getGapPos1() || 
                        val == userParameters->getGapPos2())
                    {
                        seq[i - firstRes] = '-';
                    }
                    else
                    {
                        seq[i - firstRes] = val;
                    }
                }

                //for(i = 1; i <= lastRes; i++) //- nige
                for(i = 1; i <= length; i++) //- nige
                {
                    (*gdeOutFile) << seq[i - 1];
                    //if((i % lineLength == 0) || (i == lastRes)) //- nige
                    if((i % lineLength == 0) || (i == length)) //- nige
                    {
                        (*gdeOutFile) << "\n";
                    }
                }
            }
        }
        gdeOutFile->close();   
    }
    catch(const bad_alloc& e)
    {
        gdeOutFile->close();
        cerr << "A bad_alloc exception has occured in the gdeOut function.\n"
             << e.what() << "\n";
        exit(1);
    }
    catch(VectorOutOfRange e)
    {
        gdeOutFile->close();
        cerr << "An exception has occured in the gdeOut function.\n"
             << e.what() << "\n";
        exit(1);    
    }
    catch(...)
    {
        gdeOutFile->close();
        cerr << "An exception has occured in the gdeOut function.\n";
        exit(1);    
    }     
}

void AlignmentOutput::nbrfOut(Alignment* alignPtr, outputRegion partToOutput)
{
    char residue;
    int val;
    int i, ii;
    int j, slen;    
    int lineLength;
    rangeNum rnum;
    int firstRes = partToOutput._firstRes;
    int lastRes = partToOutput._lastRes;
    int firstSeq = partToOutput._firstSeq;
    int lastSeq = partToOutput._lastSeq;
            
    try
    {
        int _maxLength = alignPtr->getMaxAlnLength();
        vector<char> sequence;
        sequence.assign(_maxLength + 1, '0');
        const SeqArray* alignment = alignPtr->getSeqArray();
        // decide the line length for this alignment - maximum is LINELENGTH 
        lineLength = PAGEWIDTH - alignPtr->getMaxNames();
        lineLength = lineLength - lineLength % 10; // round to a multiple of 10
        if (lineLength > LINELENGTH || lineLength <= 0) // Mark, 18-7-07
        {
            lineLength = LINELENGTH;
        }
        
        // Get name prefix. DL if DNA, P1 if protein
        string namePrefix;
        bool _dnaFlag = userParameters->getDNAFlag();
        namePrefix = _dnaFlag ? ">DL;" : ">P1;";
        
        //int length = lastRes - firstRes + 1; //- nige

        for(ii = firstSeq; ii <= lastSeq; ii++) 
        {
            i = alignPtr->getOutputIndex(ii - 1);
            
            (*nbrfOutFile) << namePrefix;
            
            if (!userParameters->getSeqRange()) 
            {
                (*nbrfOutFile) << alignPtr->getName(i) << "\n"
                               << alignPtr->getTitle(i) << "\n";
            }
            else 
            {
                findRangeValues(alignPtr, &rnum, firstRes, lastRes, ii);
                (*nbrfOutFile) << nameonly(alignPtr->getName(i)) << "/" << rnum.start
                               << "-"<< rnum.end << "\n"
                               << alignPtr->getTitle(i) << "\n";
            }
            slen = 0;
            //for(j = firstRes; j < firstRes + lastRes; j++) //- nige
            for(j = firstRes; j <= lastRes; j++) //- nige
            {
                // NOTE I changed this here!!!!!
                if (j <= alignPtr->getSeqLength(i))
                {
                    //val = alignPtr.getResidue(i, j);
                    val = (*alignment)[i][j]; // NOTE june29
                }
                else
                { 
                    val = -3;
                }
                if((val == -3) || (val == 253))
                {
                    break;
                }
                else if((val < 0) || (val > userParameters->getMaxAA()))
                {
                    residue = '-';
                }
                else 
                {
                    residue = userParameters->getAminoAcidCode(val);
                }
                sequence[j - firstRes] = residue;
                slen++;
            }
            for(j = 1; j <= slen; j++) 
            {
                (*nbrfOutFile) << sequence[j - 1];
                if((j % lineLength == 0) || (j == slen))
                { 
                    (*nbrfOutFile) << "\n";
                }
            }
            (*nbrfOutFile) << "*\n";
        }    
        nbrfOutFile->close();
    }
    catch(const bad_alloc& e)
    {
        nbrfOutFile->close();
        cerr << "A bad_alloc exception has occured in the nbrfOut function.\n"
             << e.what() << "\n";
        exit(1);
    }
    catch(VectorOutOfRange e)
    {
        nbrfOutFile->close();
        cerr << "An exception has occured in the nbrfOut function.\n"
             << e.what() << "\n";
        exit(1);    
    }
    catch(...)
    {
        nbrfOutFile->close();
        cerr << "An exception has occured in the nbrfOut function.\n";
        exit(1);    
    }    
}

/* 
 * The function clustalOut is used to ouput the Alignment into a clustal format
 * file.
 */
void AlignmentOutput::clustalOut(Alignment* alignPtr, outputRegion partToOutput)
{
    int firstRes = partToOutput._firstRes;
    int lastRes = partToOutput._lastRes;
    int firstSeq = partToOutput._firstSeq;
    int lastSeq = partToOutput._lastSeq;
        
    try
    {
        vector<char> seq1;    
        vector<int> seqNo;
        vector<int> printSeqNo;    
        vector<char> ss_mask1, ss_mask2;
        const SeqArray* alignment = alignPtr->getSeqArray();
        char  temp[MAXLINE];
        char c;
        int val;
        int ii, lv1, catident1[NUMRES], catident2[NUMRES], ident, chunks;
        int i, j, k, l;
        int pos, ptr;
        int lineLength;

        rangeNum rnum;
        string tmpStr = "";
        string sequenceLine = "";


        int _numSequences = alignPtr->getNumSeqs();
        // PMcG revert to Clustalw1.83 output format for name width
        const int NAMESWIDTH=std::max(std::min(MAXNAMESTODISPLAY,alignPtr->getMaxNames()) , MINNAMESTODISPLAY); 

        seqNo.assign(_numSequences + 1, 0);
        printSeqNo.assign(_numSequences + 1, 0);
    
        // check that lastSeq <= _numSequences
        if(lastSeq > _numSequences)
        {
            lastSeq = _numSequences;
        }
    
        for (i = firstSeq; i <= lastSeq; i++)
        {
            printSeqNo[i] = seqNo[i] = 0;
            for(j = 1; j < firstRes; j++) 
            {
                //val = alignPtr.getResidue(i, j);
                val = (*alignment)[i][j]; // NOTE june29
                if((val >= 0) || (val <= userParameters->getMaxAA()))
                { 
                    seqNo[i]++;
                }
            }
        }

        seq1.assign(alignPtr->getMaxAlnLength() + 1, 0);
        // Check if we have secondary structure in file 1 and if we want to output it.
        if (userParameters->getStructPenalties1() == SECST && 
            userParameters->getUseSS1() == true) 
        {
            int lengthSeq1 = alignPtr->getSeqLength(1);
            ss_mask1.assign(lengthSeq1 + 10, 0);
            vector<char>* _secStructMask1 = alignPtr->getSecStructMask1();
        
            for (i = 0; i < lengthSeq1; i++)
            {
                ss_mask1[i] = _secStructMask1->at(i);
            }
            printSecStructMask(lengthSeq1, _secStructMask1, &ss_mask1);
          
        }
        
        // Check if we have secondary structure in file 2 and if we want to output it.
        if (userParameters->getStructPenalties2() == SECST && 
            userParameters->getUseSS2() == true &&
            userParameters->getProfile2Empty() == false)
        {
            // AW added test getProfile2Empty was meant to fix bug
            // 141, but only prevents segfault, output still faulty.
            // Has to be fixed somewhere else
            int indexProfile2FirstSeq = alignPtr->getProfile1NumSeqs() + 1;
            int lengthSeqProfile2 = alignPtr->getSeqLength(indexProfile2FirstSeq);
            ss_mask2.assign(lengthSeqProfile2 + 10, 0);
            vector<char>* _secStructMask2 = alignPtr->getSecStructMask2();
            for (i = 0; i < lengthSeqProfile2; i++)
            {
                ss_mask2[i] = _secStructMask2->at(i);
            }
            printSecStructMask(lengthSeqProfile2, _secStructMask2, &ss_mask2);
        }
    
        (*clustalOutFile) << "CLUSTAL "<< userParameters->getRevisionLevel() 
                          << " multiple sequence alignment\n\n";
    
        // decide the line length for this alignment - maximum is LINELENGTH

        //PMcG  9-2-2008 make line output same as 1.83
        //lineLength = PAGEWIDTH - alignPtr->getMaxNames();
        lineLength = PAGEWIDTH - NAMESWIDTH;
        lineLength = lineLength - lineLength % 10; // round to a multiple of 10
        if (lineLength > LINELENGTH || lineLength <= 0) // Mark 18-7-07
        { 
            lineLength = LINELENGTH;
        }
    
        int length = lastRes - firstRes + 1; //- nige

        chunks = length / lineLength; //- nige
        if(length % lineLength != 0) //- nige
        {
          ++chunks;
        }
        //printf("firstRes=%d,lastRes=%d,length=%d,chunks=%d\n",
        //       firstRes, lastRes, length, chunks);

        // This will loop through each of the blocks.
        for(lv1 = 1; lv1 <= chunks; ++lv1) 
        {
            // pos is begining of chunk, ptr is the end of the chunk to be displayed.
            pos = ((lv1 - 1) * lineLength) + 1; 
            ptr = (length < pos + lineLength - 1) ? length : pos + lineLength - 1; //- nige
   
            (*clustalOutFile) << "\n";
            int _outStructPenalty = userParameters->getOutputStructPenalties();
            
            string secStructName1 = alignPtr->getSecStructName1(); // Mark 18-7-07
            
            if((int)secStructName1.size() > MAXNAMESTODISPLAY)
            {
                //secStructName1 = secStructName1.substr(0, MAXNAMESTODISPLAY);
                secStructName1 = secStructName1.substr(0, NAMESWIDTH);
            }            
            if (_outStructPenalty == 0 || _outStructPenalty == 2) 
            {
                if (userParameters->getStructPenalties1() == SECST && 
                    userParameters->getUseSS1() == true) 
                {
                    for(i = pos; i <= ptr; ++i) 
                    {
                        // Check if we can access mask position first
                        if((int)ss_mask1.size() > i + firstRes - 2)  
                        {
                            val = ss_mask1[i + firstRes - 2]; 
                            if (val == userParameters->getGapPos1() || 
                                val == userParameters->getGapPos2())
                            {
                                temp[i - pos] = '-';
                            }
                            else
                            {
                                temp[i - pos] = val;
                            }
                        }
                        else
                        {
                            break;
                        }
                    }
                    temp[i - pos] = EOS;  
                                      
                    if(userParameters->getSeqRange())
                    {
                        //(*clustalOutFile) << "!SS_" << setw(MAXNAMESTODISPLAY +
                        (*clustalOutFile) << "!SS_" << setw(NAMESWIDTH +
                                                            clusSecStructOffset)
                                          << left << secStructName1 
                                          << "  " << temp << "\n";
                    }
                    else
                    {
                        //(*clustalOutFile) << "!SS_" << setw(MAXNAMESTODISPLAY)
                        (*clustalOutFile) << "!SS_" << setw(NAMESWIDTH)
                                          << left << secStructName1 << "  " 
                                          << temp << "\n";
                    }
                }
            }
            if (_outStructPenalty == 1 || _outStructPenalty == 2) 
            {
                if (userParameters->getStructPenalties1() != NONE && 
                    userParameters->getUseSS1() == true) 
                {
                    vector<char>* _gapPenaltyMask1 = alignPtr->getGapPenaltyMask1();       
                    for(i = pos; i <= ptr; ++i) 
                    {
                        if((int)_gapPenaltyMask1->size() > i + firstRes - 2) 
                        {
                            val = _gapPenaltyMask1->at(i + firstRes - 2); 
                            if (val == userParameters->getGapPos1() || 
                                val == userParameters->getGapPos2())
                            {
                                temp[i - pos] = '-';
                            }
                            else
                            {
                                temp[i - pos] = val;
                            }
                        }
                        else
                        {
                            break;
                        }
                    }
                    temp[i - pos] = EOS;
                    
                    //(*clustalOutFile) << "!GM_" << setw(MAXNAMESTODISPLAY) << left
                    (*clustalOutFile) << "!GM_" << setw(NAMESWIDTH) << left
                                      << secStructName1 << "  "
                                      << temp << "\n";
                }
            }
            
            
            string secStructName2 = alignPtr->getSecStructName2(); // Mark 18-7-07
            //if((int)secStructName2.size() > MAXNAMESTODISPLAY)
            if((int)secStructName2.size() > NAMESWIDTH)
            {
                //secStructName2 = secStructName2.substr(0, MAXNAMESTODISPLAY);
                secStructName2 = secStructName2.substr(0, NAMESWIDTH);
            }
                                
            if (_outStructPenalty == 0 || _outStructPenalty == 2) 
            {
                if (userParameters->getStructPenalties2() == SECST && 
                    userParameters->getUseSS2() == true) 
                {
                    for(i = pos; i <= ptr; ++i) 
                    {
                        if((int)ss_mask2.size() > i + firstRes - 2)
                        {
                            val = ss_mask2[i + firstRes - 2]; 
                            if (val == userParameters->getGapPos1() || 
                                val == userParameters->getGapPos2())
                            {
                                temp[i - pos] = '-';
                            }
                            else
                            {
                                temp[i - pos] = val;
                            }
                        }
                        else
                        {
                            break;
                        }
                    }
                    temp[i - pos] = EOS;
                    
                                        
                    if (userParameters->getSeqRange())
                    {
                        //(*clustalOutFile) << "!SS_" << setw(MAXNAMESTODISPLAY +
                        (*clustalOutFile) << "!SS_" << setw(NAMESWIDTH +
                                                            clusSecStructOffset)
                                          << left << secStructName2; 
                        (*clustalOutFile) << "  " << temp << "\n";
                    }
                    else
                    {
                        //(*clustalOutFile) << "!SS_" << setw(MAXNAMESTODISPLAY)
                        (*clustalOutFile) << "!SS_" << setw(NAMESWIDTH)
                                          << left << secStructName2 << "  " 
                                          << temp << "\n";
                    }
                }
            }
            
            if (_outStructPenalty == 1 || _outStructPenalty == 2) 
            {
                if (userParameters->getStructPenalties2() != NONE && 
                    userParameters->getUseSS2() == true) 
                {
                    vector<char>* _gapPenaltyMask2 = alignPtr->getGapPenaltyMask2();
                    for(i = pos; i <= ptr; ++i) 
                    {
                        if((int)_gapPenaltyMask2->size() > i + firstRes - 2)
                        {
                            val = _gapPenaltyMask2->at(i + firstRes - 2); 
                            if (val == userParameters->getGapPos1() || 
                                val == userParameters->getGapPos2())
                            {
                                temp[i - pos] = '-';
                            }
                            else
                            {
                                temp[i - pos] = val;
                            }
                        }
                        else
                        {
                            break;
                        }
                    }
                    temp[i - pos] = EOS;
                    
                    //(*clustalOutFile) << "!GM_" << setw(MAXNAMESTODISPLAY) << left
                    (*clustalOutFile) << "!GM_" << setw(NAMESWIDTH) << left
                                      << secStructName2 << "  "
                                      << temp << "\n";
                }
            }
        
            for(ii = firstSeq; ii <= lastSeq; ++ii) 
            {
                i = alignPtr->getOutputIndex(ii - 1); 
                printSeqNo[i] = 0;
                for(j = pos; j <= ptr; ++j) 
                {
                    if (j + firstRes - 1 <= alignPtr->getSeqLength(i))
                    {
                        //val = alignPtr.getResidue(i, j + firstRes - 1);
                        val = (*alignment)[i][j + firstRes - 1]; // NOTE june29
                    }
                    else
                    { 
                        val = -3;
                    }
                    if((val == -3) || (val == 253)) 
                    {
                        break;
                    }
                    else if((val < 0) || (val > userParameters->getMaxAA()))
                    {
                        seq1[j] = '-';
                    }
                    else 
                    {
                        seq1[j] = userParameters->getAminoAcidCode(val);
                        seqNo[i]++;
                        printSeqNo[i] = 1;
                    } 
                }
                for(; j <= ptr; ++j) 
                {
                    seq1[j]='-';
                }
            
                sequenceLine = "";

                for(int index = pos; index < ptr + 1; index++)
                {
                    sequenceLine += seq1[index];
                }
                
                string seqName = alignPtr->getName(i);
                //if((int)seqName.size() > MAXNAMESTODISPLAY)
                if((int)seqName.size() > NAMESWIDTH)
                {
                    //seqName = seqName.substr(0, MAXNAMESTODISPLAY);
                    seqName = seqName.substr(0, NAMESWIDTH);
                }
                
                if (!userParameters->getSeqRange()) 
                {
                    // NOTE I made a change here from + 5, to + 6.
                    //(*clustalOutFile) << setw(alignPtr->getMaxNames() + 6) << left
                    //                  << alignPtr->getName(i);

                    //(*clustalOutFile) << setw(MAXNAMESTODISPLAY + 6) << left
                    //                  << seqName; 
                    // PMcG changed this to revert to behaviour of clustalw1.83
                    // 
                    //(*clustalOutFile) << setw(std::max(std::min(MAXNAMESTODISPLAY,alignPtr->getMaxNames()) , MINNAMESTODISPLAY) + 6)
                    (*clustalOutFile) << setw(NAMESWIDTH + 6) << left << seqName; 


                }
                else // Put the sequence range information in! 
                {
                    findRangeValues(alignPtr, &rnum, firstRes, lastRes, ii);
                    string nameOnly = nameonly(seqName);
                    std::stringstream ss;
                    std::stringstream ss2;
                    string rangeStart;
                    string rangeEnd;
                    ss << rnum.start;
                    ss >> rangeStart;
                    ss2 << rnum.end;
                    ss2 >> rangeEnd;
                    tmpStr = nameOnly + "/" + rangeStart + "-" + rangeEnd;
                    //(*clustalOutFile) << setw(MAXNAMESTODISPLAY + clusSequenceOffset) 
                    (*clustalOutFile) << setw(NAMESWIDTH + clusSequenceOffset) 
                                      << left << tmpStr;
                }
            
                (*clustalOutFile) << sequenceLine;
            
                if (userParameters->getClSeqNumbers() && printSeqNo[i])
                {
                    (*clustalOutFile) << " " << seqNo[i];
                }
            
                (*clustalOutFile) << "\n";
            }
        
            // Now print out the conservation information!
            for(i = pos; i <= ptr; ++i) 
            {
                seq1[i] = ' ';
                ident = 0;

                for(j = 1; j <= (int)strongGroup.size(); j++)
                {
                    catident1[j - 1] = 0;
                }
                for(j = 1; j <= (int)weakGroup.size(); j++) 
                {
                    catident2[j - 1] = 0;
                }
                
                for(j = firstSeq; j <= lastSeq; ++j) 
                {
                    if((i + firstRes - 1 <= alignPtr->getSeqLength(j)) &&
                       (i + firstRes - 1 <= alignPtr->getSeqLength(firstSeq)))
                    {// NOTE june29
                        if(((*alignment)[firstSeq][i + firstRes - 1] >= 0) && 
                            ((*alignment)[firstSeq][i + firstRes - 1] <=
                            userParameters->getMaxAA())) 
                        {
                            // Count how many are identical to the first sequence
                            if((*alignment)[firstSeq][i + firstRes - 1] == 
                                (*alignment)[j][i + firstRes - 1])
                            {
                                ++ident;
                            }
                            // Count how many are in the same category.
                            for(k = 1; k <= (int)strongGroup.size(); k++) 
                            {
                                for(l = 0; (c = strongGroup[k - 1][l]); l++) 
                                {
                                    int resCode = (*alignment)[j][i + firstRes - 1];
                                    if(resCode <= userParameters->getMaxAA() + 1)
                                    {
                                        if (userParameters->getAminoAcidCode(resCode) == c)
                                        {
                                            catident1[k - 1]++;
                                            break;
                                        }
                                    }
                                }
                            }
                            for(k = 1; k <= (int)weakGroup.size(); k++) 
                            {
                                for(l = 0; (c = weakGroup[k - 1][l]); l++) 
                                {//NOTE june29
                                    int resCode = (*alignment)[j][i + firstRes - 1];
                                    if(resCode <= userParameters->getMaxAA() + 1)
                                    {
                                        if (userParameters->getAminoAcidCode(resCode) == c)
                                        {
                                            catident2[k - 1]++;
                                            break;
                                        }
                                    }
                                }
                            }
                        }
                    }    
                }
            
                // Now do the conservation part for each block. 
                if(ident == lastSeq - firstSeq + 1)
                {
                    seq1[i] = '*'; // All residues the same!
                }
                else if (!userParameters->getDNAFlag()) 
                {
                    for(k = 1; k <= (int)strongGroup.size(); k++) 
                    {
                        if (catident1[k - 1] == lastSeq - firstSeq + 1) 
                        {
                            seq1[i] = ':'; // All residues member of the same category
                            break;
                        }
                    }
                    if(seq1[i] == ' ')
                    {
                        for(k = 1; k <= (int)weakGroup.size(); k++) 
                        {
                            if (catident2[k - 1] == lastSeq - firstSeq + 1) 
                            {
                                seq1[i] = '.'; // All residues member of the same category
                                break;
                            }
                        }
                    }
                }
            }

            sequenceLine = "";
            for(int index = pos; index < ptr + 1; index++)
            {
                sequenceLine += seq1[index];
            }

            //for(k = 0; k < MAXNAMESTODISPLAY + 6; k++)
            for(k = 0; k < NAMESWIDTH + 6; k++)
            { 
                (*clustalOutFile) << " ";
            }
            if(userParameters->getSeqRange()) 
            {

                (*clustalOutFile) << "          ";
            }

            (*clustalOutFile) << sequenceLine << "\n";
        }
        clustalOutFile->close();
    }
    catch(const bad_alloc& e)
    {
        clustalOutFile->close();
        cerr << "A bad_alloc exception has occured in the clustalOut function.\n"
             << e.what() << "\n";
        exit(1);
    }
    catch(VectorOutOfRange e)
    {
        clustalOutFile->close();
        cerr << "An exception has occured in the clustalOut function.\n"
             << e.what() << "\n";
        exit(1);    
    }
    catch(exception e)
    {
        clustalOutFile->close();
        cerr << "An exception has occured in the clustalOut function.\n"
             << e.what() << "\n";
        exit(1);    
    }
    catch(...)
    {
        clustalOutFile->close();
        cerr << "An exception has occured in the clustalOut function.\n";
        exit(1);    
    }    
}

/*
 * The function printSecStructMask takes in the 2 mask vectors. It makes changes to 
 * mask and puts the result in structMask. This is the one that will be printed out.
 * NOTE I have had to use (*structMask)[i] syntax. This is not good. But there is no
 * way to do random access setting in vectors otherwise.
 */
void AlignmentOutput::printSecStructMask(int prfLength, vector<char>* mask, 
                                         vector<char>* structMask)
{
    int i, j;
//    calculate the gap penalty mask from the secondary structures

    i = 0;
    try
    {
        while (i < prfLength) 
        {
            if (tolower(mask->at(i)) == 'a' || mask->at(i) == '$') 
            {
                for (j = 0; j < userParameters->getHelixEndMinus(); j++) 
                {
                    if (i + j >= prfLength || (tolower(mask->at(i+j)) != 'a'
                                    && mask->at(i+j) != '$')) 
                    {    
                        break;
                    }
                    (*structMask)[i+j] = 'a';
                }
                i += j;
                while (tolower(mask->at(i)) == 'a' || mask->at(i) == '$') 
                {
                    if (i >= prfLength) 
                        break;
                    if (mask->at(i) == '$') 
                    {
                        (*structMask)[i] = 'A';
                        i++;
                        break;
                    }
                    else 
                        (*structMask)[i] = (*mask)[i];
                    i++;
                }
                for (j = 0; j < userParameters->getHelixEndMinus(); j++) 
                {
                    if ((i-j-1>=0) && (tolower(mask->at(i-j-1)) == 'a'
                                        || mask->at(i-j-1) == '$'))
                    {
                        (*structMask)[i - j - 1] = 'a';
                    }
                }
            }
            else if (tolower(mask->at(i)) == 'b' || mask->at(i) == '%') 
            {
                for (j = 0; j < userParameters->getStrandEndMinus(); j++) 
                {
                    if (i + j >= prfLength || (tolower(mask->at(i+j)) != 'b'
                                        && mask->at(i+j) != '%')) 
                    {
                        break;
                    }
                    (*structMask)[i+j] = 'b';
                }
                i += j;
                while (tolower(mask->at(i)) == 'b' || mask->at(i) == '%') 
                {
                    if (i >= prfLength) 
                        break;
                    if (mask->at(i) == '%') 
                    {
                        (*structMask)[i] = 'B';
                        i++;
                        break;
                    }
                    else
                    { 
                        (*structMask)[i] = mask->at(i);
                    }
                    i++;
                }
                for (j = 0; j < userParameters->getStrandEndMinus(); j++) 
                {
                    if ((i - j - 1 >= 0) && (tolower(mask->at(i - j - 1)) == 'b' || 
                        mask->at(i-j-1) == '%'))
                    {
                        (*structMask)[i-j-1] = 'b';
                    }
                }
            }
            else
            { 
                i++;
            }
        }
    }
    catch(const exception& e)
    {
        // Catch all the exceptions
        cerr << "There has been an exception in the function printSecStructMask\n"
             << e.what() << "\n";
        exit(1);
    }
}

/*
 * The function findRangeValues is used to find the range to be printed out.
 *
 *
 */
void AlignmentOutput::findRangeValues(Alignment* alignPtr, rangeNum *rnum, int firstRes, 
                                int len, int firstSeq)
{
    assert(rnum);
    assert(firstRes > 0);
    assert(len > 0);
    assert(firstSeq > 0);
    
    try
    {
        int val;
        int i, ii;
        int j, slen;    

        char tmpName[FILENAMELEN + 15];
        int iStart = 0;
        int iEnd = 0; // to print sequence start-end with names
        int found = 0;
        int ngaps = 0;
        int tmpStart = 0; 
        int tmpEnd = 0;
        int ntermgaps = 0;
        int pregaps = 0;
        int tmpk = 0;
        int isRange = 0;
        int formula = 0;
        const SeqArray* alignment = alignPtr->getSeqArray();
        
        tmpName[0] = '\0';
        slen = 0;

        ii = firstSeq ;
        i = alignPtr->getOutputIndex(ii - 1); // NOTE Same as elsewhere! 
        string name = alignPtr->getName(i);
        
        if( (sscanf(name.c_str(),"%[^/]/%d-%d",tmpName, &tmpStart, &tmpEnd) == 3)) 
        {
            isRange = 1;
        }
        
        for(tmpk = 1; tmpk < firstRes; tmpk++)
        { 
            // do this irrespective of above sscanf
            //val = alignPtr.getResidue(i, tmpk);
            val = (*alignment)[i][tmpk]; // NOTE june29
            if ((val < 0) || (val > userParameters->getMaxAA()))
            { 
                //it is gap
                pregaps++;
            }
        }
        
        for(j = firstRes; j < firstRes + len; j++) 
        {
            if(j > alignPtr->getSeqLength(i))
            {
                val = -3; // Cant get it so break out.
            }
            else
            {
                //val = alignPtr.getResidue(i, j);
                val = (*alignment)[i][j]; // NOTE june29
            }
            if((val == -3) || (val == 253))
            {
                break;
            }
            else if((val < 0) || (val > userParameters->getMaxAA())) 
            {
                ngaps++;
            }
            else 
            {
                found = j;
            }
            if ( found && (iStart == 0) ) 
            {
                iStart = found;
                ntermgaps = ngaps;
            }
            slen++;
        }
        if(userParameters->getSeqRange()) 
        {
            cout << "Name : " << alignPtr->getName(i) << " "
                 << "\n  firstRes = "<< firstRes << " "
                 << "   len = " << len << " "
                 << "\n  iStart = " << iStart << " "
                 << "\n  tmpStart = " << tmpStart << " "
                 << "\n  ngaps = " << ngaps << " "
                 << "\n  pregaps = " << pregaps << " ";
            if (!isRange)
            {
                formula = iStart - pregaps;
            }
            else
            {
                formula = iStart - pregaps +  ( tmpStart == 1 ? 0: tmpStart-1) ;
            }

            cout << "\n\nsuggestion  iStart - pregaps + tmpStart - ntermgaps = "
                 << iStart << " - " << pregaps << " + " << tmpStart << " - " << ntermgaps
                 << " formula " << formula << " ";
        }
        else 
        {
            cerr << "\n no range found .... strange,  iStart = " << iStart;
            formula = 1;
        }
        if (pregaps == firstRes - 1) // all gaps -  now the conditions........ 
        {
            formula = tmpStart ; // keep the previous start... 
        }
        formula = (formula <= 0) ? 1: formula;
        if (pregaps == 0 && tmpStart == 0) 
        {
            formula = firstRes;
        }
        iEnd = formula + len - ngaps -1;
        rnum->start = formula;
        rnum->end = iEnd;
        cout << "\n check... " << alignPtr->getName(i) << " " << rnum->start 
             << " - " << rnum->end;
        cout << " Done checking.........";
    }
    catch(VectorOutOfRange e)
    {
        cerr << "An exception has occured in the findRangeValues function.\n"
             << e.what() << "\n";
        exit(1);    
    }    
    catch(...)
    {
        cerr << "An exception has occured in findRangeValues function\n";
        exit(1);
    }
}

string AlignmentOutput::nameonly(string s)
{
    string tmp;
    try
    {
        int i = 0;

        while (i < (int)s.size()) 
        {
            if(s.at(i) != '/')
            {
                tmp += s.at(i);
                i++;
            }
            else
            {
                break;
            }
        }
    }
    catch(const exception& e)
    {
        cerr << "An exception has occured in the function nameonly\n"
             << e.what();
        exit(1);
    }
    return tmp;
}

int AlignmentOutput::SeqGCGCheckSum(vector<char>* seq, int length)
{
    int i;
    long check;
    int seqResIndex = 1;
    //int seqLength = seq->size();
    
    for (i = 0, check = 0; i < length; i++, seqResIndex++)
    {
        char _toUpperCase = (*seq)[seqResIndex];
        check += ((i % 57) + 1) * toupper(_toUpperCase);
    }

    return (check % 10000);
}

void AlignmentOutput::showAlign()
{
    //FILE *file;
    int  numLines;
    char temp[MAXLINE + 1];
    temp[0] = '\0';
    string fileName;
    string answer;
    
    if(userParameters->getOutputClustal()) 
    {
        fileName = clustalOutName;
    }
    else if(userParameters->getOutputNbrf()) 
    {
        fileName = nbrfOutName;
    }
    else if(userParameters->getOutputGCG()) 
    {
        fileName = gcgOutName;
    }
    else if(userParameters->getOutputPhylip()) 
    {
        fileName = phylipOutName;
    }
    else if(userParameters->getOutputGde()) 
    {
        fileName = gdeOutName;
    }
    else if(userParameters->getOutputNexus()) 
    {
        fileName = nexusOutName;
    }
    else if(userParameters->getOutputFasta()) 
    {
        fileName = fastaOutName;
    }
    else
    {
        return; // No file output type!
    }

    ifstream _fileIn;
    _fileIn.open(fileName.c_str(), ios::in);        
    _fileIn.seekg(0, std::ios::beg); // start at the beginning
    
    cout << "\n\n";
    numLines = 0;

    while(_fileIn.getline(temp, MAXLINE + 1, '\n')) 
    {
       //fputs(temp,stdout);
       cout << temp << "\n";
       ++numLines;
       if(numLines >= PAGE_LEN) 
       {
           cout << "\n";
           utilityObject->getStr(string("Press [RETURN] to continue or  X  to stop"), answer);
           if(toupper(answer[0]) == 'X') 
           {
               _fileIn.close();
               return;
           }
           else
           {
               numLines = 0;
           }
       }
    }
    _fileIn.close();
    cout << "\n";
    utilityObject->getStr(string("Press [RETURN] to continue"), answer);
}

}