/* massXpert - the true massist's program.
   --------------------------------------
   Copyright(C) 2006,2007 Filippo Rusconi

   http://www.massxpert.org/massXpert

   This file is part of the massXpert project.

   The massxpert project is the successor to the "GNU polyxmass"
   project that is an official GNU project package(see
   www.gnu.org). The massXpert project is not endorsed by the GNU
   project, although it is released ---in its entirety--- under the
   GNU General Public License. A huge part of the code in massXpert
   is actually a C++ rewrite of code in GNU polyxmass. As such
   massXpert was started at the Centre National de la Recherche
   Scientifique(FRANCE), that granted me the formal authorization to
   publish it under this Free Software License.

   This software is free software; you can redistribute it and/or
   modify it under the terms of the GNU  General Public
   License version 3, as published by the Free Software Foundation.
   

   This software is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
   General Public License for more details.

   You should have received a copy of the GNU General Public License
   along with this software; if not, write to the

   Free Software Foundation, Inc.,

   51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
*/


/////////////////////// Local includes
#include "fragmenter.hpp"


namespace massXpert
{

  Fragmenter::Fragmenter(Polymer *polymer,
			  const PolChemDef *polChemDef, 
			  const QList<FragOptions *> &fragOptionList,
			  const CalcOptions &calcOptions, 
			  const IonizeRule &ionizeRule)
    : mp_polymer(polymer), 
      mp_polChemDef(polChemDef),
      m_calcOptions(calcOptions),
      m_ionizeRule(ionizeRule)
  {
    Q_ASSERT(mp_polymer && mp_polChemDef);
  
    for (int iter = 0; iter < fragOptionList.size(); ++iter)
      {
	FragOptions *fragOptions = fragOptionList.at(iter);
	FragOptions *newFragOptions = new FragOptions(*fragOptions);
      
	m_fragOptionList.append(newFragOptions);
      }
  }


  Fragmenter::Fragmenter(const Fragmenter &other)
    : mp_polymer(other.mp_polymer), 
      mp_polChemDef(other.mp_polChemDef), 
      m_calcOptions(other.m_calcOptions),
      m_ionizeRule(other.m_ionizeRule)
  {
    Q_ASSERT(mp_polymer && mp_polChemDef);
  
    for (int iter = 0; iter < other.m_fragOptionList.size(); ++iter)
      {
	FragOptions *fragOptions = other.m_fragOptionList.at(iter);
	FragOptions *newFragOptions = new FragOptions(*fragOptions);
      
	m_fragOptionList.append(newFragOptions);
      }
  }


  Fragmenter::~Fragmenter()
  {
    freeOligomerLists();
    qDeleteAll(m_fragOptionList.begin(), m_fragOptionList.end());
    m_fragOptionList.clear();
  }


  // Takes ownership of the parameter.
  void
  Fragmenter::addFragOptions(FragOptions *fragOptions)
  {
    Q_ASSERT(fragOptions);
  
    m_fragOptionList.append(fragOptions);
  }


  QList<OligomerList *> &
  Fragmenter::oligomerLists()
  {
    return m_oligomerLists;
  }


  bool
  Fragmenter::fragment()
  {
    // If the polymer sequence is empty, just return.
    if (!mp_polymer->size())
      return true;
  
    // Ensure that the list of fragmentation options is not empty.

    if (!m_fragOptionList.size())
      {
	qDebug() << __FILE__ << __LINE__
		  << "List of fragmentation options is empty !";
      
	return false;
      }
  
    //   qDebug() << __FILE__ << __LINE__
    //  	    << "number of fragmentation specifications:" 
    //  	    << m_fragOptionList.size();
  
    // For each fragmentation options instance in the list, perform the
    // required fragmentation.

    for (int iter = 0; iter < m_fragOptionList.size(); ++iter)
      {
	FragOptions *fragOptions = m_fragOptionList.at(iter);
      
	if(fragOptions->fragEnd() == MXT_FRAG_END_NONE)
	  {
	    if (fragmentEndNone(*fragOptions) == -1)
	      return false;
	  }
	else if (fragOptions->fragEnd() == MXT_FRAG_END_LEFT)
	  {
	    if (fragmentEndLeft(*fragOptions) == -1)
	      return false;
	  }
	else if (fragOptions->fragEnd() == MXT_FRAG_END_RIGHT)
	  {
	    if (fragmentEndRight(*fragOptions) == -1)
	      return false;
	  }
	else
	  Q_ASSERT(0); 
      }

    return true;
  }


  int
  Fragmenter::fragmentEndNone(FragOptions &fragOptions) 
  {
    int count = 0;
  
    OligomerList *oligomerList;
    oligomerList = new OligomerList(fragOptions.name(), mp_polymer);

    // We are generating fragments that are made of a single monomer,
    // like in the proteinaceous world we have the immonium ions.

    CalcOptions localOptions = m_calcOptions;
  
    for (int iter = fragOptions.startIndex(); 
	 iter < fragOptions.endIndex() + 1; ++iter)
      {
	// We create an oligomer which is not ionized(false) but that
	// bears the default ionization rule, because this oligomer
	// might be later used in places where the ionization rule has
	// to be valid. For example, one drag and drop operation might
	// copy this oligomer into a mzLab dialog window where its
	// ionization rule validity might be challenged. Because this
	// fragmentation oligomer will bear only its intrinsic 1
	// charge, we should set the level member of the ionization to
	// 0.
	
	IonizeRule ionizeRule(m_ionizeRule);
	ionizeRule.setLevel(0);
	
	FragmentOligomer *oligomer1 = 
	  new FragmentOligomer(mp_polymer, "NOT_SET", fragOptions.name(),
				Ponderable(), 
				ionizeRule, false,
				iter, iter);
	
	localOptions.setCapping(MXT_CAP_NONE);
      

	// The calculation of the masses is not performed by applying
	// any ionization rule, as the fragmentation rule itself will
	// yield an ionized species(it should be singly-charged). Thus
	// we pass an invalid ionizeRule object(upon creation, an
	// ionize rule is invalid).
	IonizeRule rule;

	if(!oligomer1->calculateMasses(&localOptions, &rule))
	  {
	    delete oligomer1;

	    return -1;
	  }

	const QList<Atom *> &refList = mp_polChemDef->atomList();
      
	// It is at this precise moment that the generated oligomer will
	// have a single charge onto it: by the accounting of the
	// fragmentation rule's formula.

	if(!fragOptions.formula().isEmpty())    
	  if (!fragOptions.Formula::accountMasses(refList, oligomer1))
	    {
	      delete oligomer1;

	      return -1;
	    }
      
	// At this moment, the new fragment might be challenged for
	// the fragmented monomer's contribution. For example, in
	// nucleic acids, it happens that during a fragmentation, the
	// base of the fragmented monomer is decomposed and goes
	// away. This is implemented in massXpert with the ability to
	// tell the fragmenter that upon fragmentation the mass of the
	// monomer is to be removed. The skeleton mass is then added
	// to the formula of the fragmentation pattern.

	int monomerContrib = fragOptions.monomerContribution();

	if(monomerContrib)
	  {
	    const Monomer *monomer = mp_polymer->at(iter);
	  
	    QString formula = monomer->formula();
	  		
	    if (!monomer->accountMasses(&oligomer1->rmono(),
					 &oligomer1->ravg(),
					 monomerContrib))
	      {
		delete oligomer1;
	      
		return -1;
	      }
	  }
      
	// At this point we should check if the fragmentation
	// specification includes fragmentation rules that apply to this
	// fragment.

	for(int jter = 0; jter < fragOptions.ruleList().size(); ++jter)
	  {
	    // The accounting of the fragrule is performed on a
	    // singly-charged oligomer, as defined by the fragmentation
	    // formula. Later, we'll have to take into account the fact
	    // that the user might want to calculate fragment m/z with
	    // z>1.

	    FragRule *fragRule = fragOptions.ruleList().at(jter);
	  
	    if (!accountFragRule(fragRule, true, iter, MXT_FRAG_END_NONE, 0))
	      continue;
	  
	    // Each fragrule triggers the creation of a new oligomer.

	    FragmentOligomer *oligomer2 = 
	      new FragmentOligomer(*oligomer1);

	    if (!accountFragRule(fragRule, false, iter, 
				  MXT_FRAG_END_NONE, oligomer2))
	      {
		delete oligomer1;
		delete oligomer2;

		return -1;
	      }

	    // At this point we have the fragment oligomer that bears
	    // the "default" fragmentation specification
	    // formula-driven charge, that is z=1.
	    
	    FragmentOligomer *newOligomer = 
	      new FragmentOligomer(*oligomer2);
	  
	    QString name = QString("%1#%2#(%3)#z=1")
	      .arg(fragOptions.name())
	      .arg(mp_polymer->at(iter)->code())
	      .arg(fragRule->name());
	  
	    newOligomer->setName(name);
	  
	    oligomerList->append(newOligomer);
	  
	    ++count;
	  
	    // However, do not forget that the user might ask for
	    // fragments that bear more than the single charge that was
	    // intrinsically computed within the formula of the
	    // fragmentation specification. Thus create as many new
	    // oligomers as needed for the different charge levels asked
	    // by the user.  Because the ionization changes the values
	    // in the oligomer, and we need a new oligomer each time, we
	    // duplicate the oligomer each time we need it.

	    int startIonizeLevel = fragOptions.startIonizeLevel();
	    int endIonizeLevel = fragOptions.endIonizeLevel();
      
	    for (int kter = startIonizeLevel; kter < endIonizeLevel; ++kter)
	      {
		IonizeRule ionizeRule(m_ionizeRule);
		ionizeRule.setLevel(kter);

		newOligomer = new FragmentOligomer(*oligomer2);
	  
		// If the result of the call below is -1, then that
		// means that there was an error and we should return
		// immediately. If it is 0, then that means that no
		// error was encountered, but that no actual ionization
		// took place, so we need not take into account the
		// oligomer.

		int res = newOligomer->ionize(ionizeRule);
	      
		if(res == -1)
		  {
		    delete newOligomer;
		    delete oligomer2;
		    delete oligomer1;
		  
		    return -1;
		  }
		else if (res == 0)
		  {
		    delete newOligomer;
		
		    continue;
		  }
	      
		// At this point the ionization did indeed perform
		// something interesting, continue constructing the
		// resulting oligomer.

		QString name = QString("%1#%2#(%3)#z=%4")
		  .arg(fragOptions.name())
		  .arg(mp_polymer->at(iter)->code())
		  .arg(fragRule->name())
		  .arg(newOligomer->charge());

		newOligomer->setName(name);
	      
		oligomerList->append(newOligomer);

		++count;
	      }
	  
	    // We can delete the template oligomer.
	    delete oligomer2;
	  }
	// End of 
	// for (int jter = 0; jter < fragOptions.ruleList().size(); ++jter)

	// We are here because of two reasons:
  
	// 1. because the array of fragRules did not contain any
	// 1. fragRule, in which case we still have to validate and
	// 1. terminate the oligomer1;
  
	// 2. because we finished dealing with fragRules, in which case
	// 2. we do not continue the creation of oligomer1, because if
	// 2. there were fragrules(even if none matched) the oligomer_1
	// 2. is not considered valid.
  
	if(!fragOptions.ruleList().size())
	  {
	    // At this point we have the fragment oligomer. However, do
	    // not forget that the user might ask for fragments that
	    // bear more than the single charge that was intrinsically
	    // computed within the formula of the fragmentation
	    // specification. 

	    // So, first create an oligomer with the "default"
	    // fragmentation specification-driven 1+ charge.  

	    FragmentOligomer *newOligomer = 
	      new FragmentOligomer(*oligomer1);
	  
	    QString name = QString("%1#%2#z=1")
	      .arg(fragOptions.name())
	      .arg(mp_polymer->at(iter)->code());
	  
	    newOligomer->setName(name);
	  
	    oligomerList->append(newOligomer);
	  
	    ++count;
	  
	    // And then create as many new oligomers as needed for the
	    // different charge levels asked by the user.  Because the
	    // ionization changes the values in the oligomer, and we
	    // need a new oligomer each time, we duplicate the oligomer
	    // each time we need it.

	    int startIonizeLevel = fragOptions.startIonizeLevel();
	    int endIonizeLevel = fragOptions.endIonizeLevel();

	    for (int kter = startIonizeLevel; kter < endIonizeLevel; ++kter)
	      {
		IonizeRule ionizeRule(m_ionizeRule);
		ionizeRule.setLevel(kter);

		FragmentOligomer *newOligomer = 
		  new FragmentOligomer(*oligomer1);
	  
		// If the result of the call below is -1, then that
		// means that there was an error and we should return
		// immediately. If it is 0, then that means that no
		// error was encountered, but that no actual ionization
		// took place, so we need not take into account the
		// oligomer.

		int res = newOligomer->ionize(ionizeRule);
	      
		if(res == -1)
		  {
		    delete newOligomer;
		    delete oligomer1;
		  
		    return -1;
		  }
		else if (res == 0)
		  {
		    delete newOligomer;
		
		    continue;
		  }
	      
		// At this point the ionization did indeed perform
		// something interesting, continue constructing the
		// resulting oligomer.

		QString name = QString("%1#%2#z=%3")
		  .arg(fragOptions.name())
		  .arg(mp_polymer->at(iter)->code())
		  .arg(newOligomer->charge());
	  
		newOligomer->setName(name);
	  
		oligomerList->append(newOligomer);
	  
		++count;
	      }

	    // We can delete the template oligomer.
	    delete oligomer1;
	  }
	else
	  delete oligomer1;
      }
    // End of 
    //   for (int iter = fragOptions.startIndex(); 
    //   fragOptions.endIndex() + 1; ++iter)
  
    m_oligomerLists.append(oligomerList);
  
    return count;
  }


  int
  Fragmenter::fragmentEndLeft(FragOptions &fragOptions)
  {
    int count = 0;
    int number = 0;
  
    OligomerList *oligomerList;
    oligomerList = new OligomerList(fragOptions.name(), mp_polymer);

    static Ponderable ponderable;
    ponderable.clearMasses();

    const QList<Atom *> &refList = mp_polChemDef->atomList();

    for (int iter = fragOptions.startIndex(); 
	 iter < fragOptions.endIndex(); ++iter, ++number)
      {
	const Monomer *monomer = mp_polymer->at(iter);
      
	monomer->accountMasses(&ponderable, 1);
      
	Ponderable ponderableTemp(ponderable);
      
	if(!fragOptions.formula().isEmpty())    
	  if (!fragOptions.Formula::accountMasses(refList, &ponderableTemp))
	    {
	      return -1;
	    }
      
	Formula formula = mp_polChemDef->leftCap();
      
	formula.accountMasses(refList, &ponderableTemp, 1);
      
	if(m_calcOptions.polymerEntities() & 
	    MXT_POLYMER_CHEMENT_LEFT_END_MODIF && !fragOptions.startIndex())
	  Polymer::accountEndModifMasses 
	   (mp_polymer,
	     MXT_POLYMER_CHEMENT_LEFT_END_MODIF,
	     &ponderableTemp);
      
	// It is at this precise moment that the generated oligomer will
	// have a single charge onto it: by the accounting of the
	// fragmentation rule's formula as computed above in the
	// ponderableTemp instance used to construct the oligomer.

	// We create an oligomer which is not ionized(false) but that
	// bears the default ionization rule, because this oligomer
	// might be later used in places where the ionization rule has
	// to be valid. For example, one drag and drop operation might
	// copy this oligomer into a mzLab dialog window where its
	// ionization rule validity might be challenged. Because this
	// fragmentation oligomer will bear only its intrinsic 1
	// charge, we should set the level member of the ionization to
	// 0.
	
	IonizeRule ionizeRule(m_ionizeRule);
	ionizeRule.setLevel(0);
	
	FragmentOligomer *oligomer1 = 
	  new FragmentOligomer(mp_polymer, "NOT_SET", fragOptions.name(),
				ponderableTemp, 
				ionizeRule, false,
				fragOptions.startIndex(), iter);
	
	// At this moment, the new fragment might be challenged for the
	// fragmented monomer's side chain contribution. For example, in
	// nucleic acids, it happens that during a fragmentation, the
	// base of the fragmented monomer is decomposed and goes
	// away. This is implemented in massXpert with the ability to
	// tell the fragmenter that upon fragmentation the mass of the
	// monomer is to be removed. The skeleton mass is then added to
	// the formula of the fragmentation pattern.

	int monomerContrib = fragOptions.monomerContribution();
      
	if(monomerContrib)
	  {
	    const Monomer *monomer = mp_polymer->at(iter);
	  
	    QString formula = monomer->formula();
	  		
	    if (!monomer->accountMasses(&oligomer1->rmono(),
					 &oligomer1->ravg(),
					 monomerContrib))
	      {
		delete oligomer1;
	      
		return -1;
	      }
	  }
      
      
	// At this point we should check if the fragmentation
	// specification includes fragmentation rules that apply to this
	// fragment.

	bool success = false;

	for(int jter = 0; jter < fragOptions.ruleList().size(); ++jter)
	  {
	    // The accounting of the fragrule is performed on a
	    // singly-charged oligomer, as defined by the fragmentation
	    // formula. Later, we'll have to take into account the fact
	    // that the user might want to calculate fragment m/z with
	    // z>1.

	    FragRule *fragRule = fragOptions.ruleList().at(jter);
	  
	    if (!accountFragRule(fragRule, true, iter, MXT_FRAG_END_LEFT, 0))
	      continue;
	  
	    // Each fragrule triggers the creation of a new oligomer.

	    FragmentOligomer *oligomer2 = 
	      new FragmentOligomer(*oligomer1);

	    if (!accountFragRule(fragRule, false, iter, 
				  MXT_FRAG_END_LEFT, oligomer2))
	      {
		delete oligomer1;
		delete oligomer2;

		return -1;
	      }

	    // At this point we have the fragment oligomer that bears
	    // the "default" fragmentation specification formula-driven
	    // charge. Let's create one such fragment oligomer for the
	    // record(see the 1 charge for the oligomer).

	    FragmentOligomer *newOligomer = 
	      new FragmentOligomer(*oligomer2);
	  
	    QString name = QString("%1#%2#(%3)#z=1")
	      .arg(fragOptions.name())
	      .arg(number + 1)
	      .arg(fragRule->name());
	  
	    newOligomer->setName(name);
	  
	    oligomerList->append(newOligomer);
	  
	    ++count;
	  
	    success = true;
	  
	    // However, do not forget that the user might ask for
	    // fragments that bear more than the single charge that was
	    // intrinsically computed within the formula of the
	    // fragmentation specification. Thus create as many new
	    // oligomers as needed for the different charge levels asked
	    // by the user.  Because the ionization changes the values
	    // in the oligomer, and we need a new oligomer each time, we
	    // duplicate the oligomer each time we need it.

	    int startIonizeLevel = fragOptions.startIonizeLevel();
	    int endIonizeLevel = fragOptions.endIonizeLevel();
      
	    for (int kter = startIonizeLevel; kter < endIonizeLevel; ++kter)
	      {
		IonizeRule ionizeRule(m_ionizeRule);
		ionizeRule.setLevel(kter);

		newOligomer = new FragmentOligomer(*oligomer2);
	  
		// If the result of the call below is -1, then that
		// means that there was an error and we should return
		// immediately. If it is 0, then that means that no
		// error was encountered, but that no actual ionization
		// took place, so we need not take into account the
		// oligomer.

		int res = newOligomer->ionize(ionizeRule);
	      
		if(res == -1)
		  {
		    delete newOligomer;
		    delete oligomer2;
		    delete oligomer1;
		  
		    return -1;
		  }
		else if (res == 0)
		  {
		    delete newOligomer;
		
		    continue;
		  }
	      
		// At this point the ionization did indeed perform
		// something interesting, continue constructing the
		// resulting oligomer.

		QString name = QString("%1#%2#(%3)#z=%4")
		  .arg(fragOptions.name())
		  .arg(number + 1)
		  .arg(fragRule->name())
		  .arg(newOligomer->charge());

		newOligomer->setName(name);
	      
		oligomerList->append(newOligomer);

		++count;
	      }
	  
	    // We can delete the template oligomer.
	    delete oligomer2;
	  }
	// End of 
	// for (int jter = 0; jter < fragOptions.ruleList().size(); ++jter)

	// We are here because of two reasons:

	// 1. because the array of fragRules did not contain any
	// 1. fragRule, in which case we still have to validate and
	// 1. terminate the oligomer1;

	// 2. because we finished dealing with fragRules, in which case
	// 2. we add oligomer1 to the list of fragments, only if no one
	// 2. of the fragrules analyzed above gave a successfully
	// 2. generated fragment(success is false).

	if(!fragOptions.ruleList().size())
	  {
	    // At this point we have the fragment oligomer. However, do
	    // not forget that the user might ask for fragments that
	    // bear more than the single charge that was intrinsically
	    // computed within the formula of the fragmentation
	    // specification. 

	    // So, first create an oligomer with the "default"
	    // fragmentation specification-driven 1+ charge.  

	    FragmentOligomer *newOligomer = 
	      new FragmentOligomer(*oligomer1);
	  
	    QString name = QString("%1#%2#z=1")
	      .arg(fragOptions.name())
	      .arg(number + 1);
	  
	    newOligomer->setName(name);
	  
	    oligomerList->append(newOligomer);
	  
	    ++count;
	  
	    // And then create as many new oligomers as needed for the
	    // different charge levels asked by the user.  Because the
	    // ionization changes the values in the oligomer, and we
	    // need a new oligomer each time, we duplicate the oligomer
	    // each time we need it.

	    int startIonizeLevel = fragOptions.startIonizeLevel();
	    int endIonizeLevel = fragOptions.endIonizeLevel();

	    for (int kter = startIonizeLevel; kter < endIonizeLevel; ++kter)
	      {
		IonizeRule ionizeRule(m_ionizeRule);
		ionizeRule.setLevel(kter);

		FragmentOligomer *newOligomer = 
		  new FragmentOligomer(*oligomer1);
	  
		// If the result of the call below is -1, then that
		// means that there was an error and we should return
		// immediately. If it is 0, then that means that no
		// error was encountered, but that no actual ionization
		// took place, so we need not take into account the
		// oligomer.

		int res = newOligomer->ionize(ionizeRule);
	      
		if(res == -1)
		  {
		    delete newOligomer;
		    delete oligomer1;
		  
		    return -1;
		  }
		else if (res == 0)
		  {
		    delete newOligomer;
		
		    continue;
		  }
	      
		// At this point the ionization did indeed perform
		// something interesting, continue constructing the
		// resulting oligomer.

		QString name = QString("%1#%2#z=%3")
		  .arg(fragOptions.name())
		  .arg(number + 1)
		  .arg(newOligomer->charge());
	  
		newOligomer->setName(name);
	  
		oligomerList->append(newOligomer);
	  
		++count;
	      }

	    // We can delete the template oligomer.
	    delete oligomer1;
	  }
	else
	  {
	    if (success)
	      delete oligomer1;
	    else
	      {
		// At this point we have the fragment oligomer. However, do
		// not forget that the user might ask for fragments that
		// bear more than the single charge that was intrinsically
		// computed within the formula of the fragmentation
		// specification. 

		// So, first create an oligomer with the "default"
		// fragmentation specification-driven 1+ charge.  

		FragmentOligomer *newOligomer = 
		  new FragmentOligomer(*oligomer1);
	  
		QString name = QString("%1#%2#z=1")
		  .arg(fragOptions.name())
		  .arg(number + 1);
	  
		newOligomer->setName(name);
	  
		oligomerList->append(newOligomer);
	  
		++count;
	  
		// And then create as many new oligomers as needed for the
		// different charge levels asked by the user.  Because the
		// ionization changes the values in the oligomer, and we
		// need a new oligomer each time, we duplicate the oligomer
		// each time we need it.

		int startIonizeLevel = fragOptions.startIonizeLevel();
		int endIonizeLevel = fragOptions.endIonizeLevel();

		for(int kter = startIonizeLevel; kter < endIonizeLevel; ++kter)
		  {
		    IonizeRule ionizeRule(m_ionizeRule);
		    ionizeRule.setLevel(kter);

		    FragmentOligomer *newOligomer = 
		      new FragmentOligomer(*oligomer1);
	  
		    // If the result of the call below is -1, then that
		    // means that there was an error and we should return
		    // immediately. If it is 0, then that means that no
		    // error was encountered, but that no actual ionization
		    // took place, so we need not take into account the
		    // oligomer.

		    int res = newOligomer->ionize(ionizeRule);
	      
		    if (res == -1)
		      {
			delete newOligomer;
			delete oligomer1;
		  
			return -1;
		      }
		    else if (res == 0)
		      {
			delete newOligomer;
		
			continue;
		      }
	      
		    // At this point the ionization did indeed perform
		    // something interesting, continue constructing the
		    // resulting oligomer.

		    QString name = QString("%1#%2#z=%3")
		      .arg(fragOptions.name())
		      .arg(number + 1)
		      .arg(newOligomer->charge());
	  
		    newOligomer->setName(name);
	  
		    oligomerList->append(newOligomer);
	  
		    ++count;
		  }
	      
		// We can delete the template oligomer.
		delete oligomer1;
	      }
	  }
      }
    // End of
    //   for (int iter = fragOptions.startIndex(); 
    //   iter < fragOptions.endIndex() + 1; ++iter, ++count)
  
    m_oligomerLists.append(oligomerList);
  
    return count;
  }
 

  int
  Fragmenter::fragmentEndRight(FragOptions &fragOptions)
  {
    int count = 0;
    int number = 0;
  
    OligomerList *oligomerList;
    oligomerList = new OligomerList(fragOptions.name(), mp_polymer);

    static Ponderable ponderable;
    ponderable.clearMasses();
  
    const QList<Atom *> &refList = mp_polChemDef->atomList();

    for (int iter = fragOptions.endIndex(); 
	 iter > fragOptions.startIndex(); --iter, ++number)
      {
	const Monomer *monomer = mp_polymer->at(iter);
      
	monomer->accountMasses(&ponderable, 1);
      
	Ponderable ponderableTemp(ponderable);
      
	if(!fragOptions.formula().isEmpty())    
	  if (!fragOptions.Formula::accountMasses(refList, &ponderableTemp))
	    {
	      return -1;
	    }
      
	Formula formula = mp_polChemDef->rightCap();
      
	formula.accountMasses(refList, &ponderableTemp, 1);
      
	if(m_calcOptions.polymerEntities() & 
	    MXT_POLYMER_CHEMENT_RIGHT_END_MODIF && 
	    fragOptions.endIndex() == mp_polymer->size() - 1)
	  Polymer::accountEndModifMasses 
	   (mp_polymer, 
	     MXT_POLYMER_CHEMENT_RIGHT_END_MODIF,
	     &ponderableTemp);
      
	// It is at this precise moment that the generated oligomer will
	// have a single charge onto it: by the accounting of the
	// fragmentation rule's formula as computed above in the
	// ponderableTemp instance used to construct the oligomer.

	// We create an oligomer which is not ionized(false) but that
	// bears the default ionization rule, because this oligomer
	// might be later used in places where the ionization rule has
	// to be valid. For example, one drag and drop operation might
	// copy this oligomer into a mzLab dialog window where its
	// ionization rule validity might be challenged. Because this
	// fragmentation oligomer will bear only its intrinsic 1
	// charge, we should set the level member of the ionization to
	// 0.
	
	IonizeRule ionizeRule(m_ionizeRule);
	ionizeRule.setLevel(0);
	
	FragmentOligomer *oligomer1 = 
	  new FragmentOligomer(mp_polymer, "NOT_SET", fragOptions.name(),
				ponderableTemp, 
				ionizeRule, false, 
				iter, fragOptions.endIndex());

	// At this moment, the new fragment might be challenged for the
	// fragmented monomer's side chain contribution. For example, in
	// nucleic acids, it happens that during a fragmentation, the
	// base of the fragmented monomer is decomposed and goes
	// away. This is implemented in massXpert with the ability to
	// tell the fragmenter that upon fragmentation the mass of the
	// monomer is to be removed. The skeleton mass is then added to
	// the formula of the fragmentation pattern.

	int monomerContrib = fragOptions.monomerContribution();
      
	if(monomerContrib)
	  {
	    const Monomer *monomer = mp_polymer->at(iter);
	  
	    QString formula = monomer->formula();
	  		
	    if (!monomer->accountMasses(&oligomer1->rmono(),
					 &oligomer1->ravg(),
					 monomerContrib))
	      {
		delete oligomer1;
	      
		return -1;
	      }
	  }

	// At this point we should check if the fragmentation
	// specification includes fragmentation rules that apply to this
	// fragment.

	bool success = false;
      
	for(int jter = 0; jter < fragOptions.ruleList().size(); ++jter)
	  {
	    // The accounting of the fragrule is performed on a
	    // singly-charged oligomer, as defined by the fragmentation
	    // formula. Later, we'll have to take into account the fact
	    // that the user might want to calculate fragment m/z with
	    // z>1.

	    FragRule *fragRule = fragOptions.ruleList().at(jter);
	  
	    if (!accountFragRule(fragRule, true, iter, MXT_FRAG_END_RIGHT, 0))
	      continue;
	  
	    // Each fragrule triggers the creation of a new oligomer.

	    FragmentOligomer *oligomer2 = 
	      new FragmentOligomer(*oligomer1);

	    if (!accountFragRule(fragRule, false, iter, 
				  MXT_FRAG_END_RIGHT, oligomer2))
	      {
		delete oligomer1;
		delete oligomer2;

		return -1;
	      }

	    // At this point we have the fragment oligomer that bears
	    // the "default" fragmentation specification formula-driven
	    // charge. Let's create one such fragment oligomer for the
	    // record(see the 1 charge for the oligomer).

	    FragmentOligomer *newOligomer = 
	      new FragmentOligomer(*oligomer2);
	  
	    QString name = QString("%1#%2#(%3)#z=1")
	      .arg(fragOptions.name())
	      .arg(number + 1)
	      .arg(fragRule->name());
	  
	    newOligomer->setName(name);
	  
	    oligomerList->append(newOligomer);
	  
	    ++count;
	  
	    success = true;
	  
	    // However, do not forget that the user might ask for
	    // fragments that bear more than the single charge that was
	    // intrinsically computed within the formula of the
	    // fragmentation specification. Thus create as many new
	    // oligomers as needed for the different charge levels asked
	    // by the user.  Because the ionization changes the values
	    // in the oligomer, and we need a new oligomer each time, we
	    // duplicate the oligomer each time we need it.

	    int startIonizeLevel = fragOptions.startIonizeLevel();
	    int endIonizeLevel = fragOptions.endIonizeLevel();
      
	    for (int kter = startIonizeLevel; kter < endIonizeLevel; ++kter)
	      {
		IonizeRule ionizeRule(m_ionizeRule);
		ionizeRule.setLevel(kter);

		newOligomer = new FragmentOligomer(*oligomer2);
	  
		// If the result of the call below is -1, then that
		// means that there was an error and we should return
		// immediately. If it is 0, then that means that no
		// error was encountered, but that no actual ionization
		// took place, so we need not take into account the
		// oligomer.

		int res = newOligomer->ionize(ionizeRule);
	      
		if(res == -1)
		  {
		    delete newOligomer;
		    delete oligomer2;
		    delete oligomer1;
		  
		    return -1;
		  }
		else if (res == 0)
		  {
		    delete newOligomer;
		
		    continue;
		  }
	      
		// At this point the ionization did indeed perform
		// something interesting, continue constructing the
		// resulting oligomer.

		QString name = QString("%1#%2#(%3)#z=%4")
		  .arg(fragOptions.name())
		  .arg(number + 1)
		  .arg(fragRule->name())
		  .arg(newOligomer->charge());

		newOligomer->setName(name);
	      
		oligomerList->append(newOligomer);

		++count;
	      }
	  
	    // We can delete the template oligomer.
	    delete oligomer2;
	  }
	// End of 
	// for (int jter = 0; jter < fragOptions.ruleList().size(); ++jter)
  
	// We are here because of two reasons:
  
	// 1. because the array of fragRules did not contain any
	// 1. fragRule, in which case we still have to validate and
	// 1. terminate the oligomer1;

	// 2. because we finished dealing with fragRules, in which case
	// 2. we add oligomer1 to the list of fragments, only if no one
	// 2. of the fragrules analyzed above gave a successfully
	// 2. generated fragment(success is false).

	if(!fragOptions.ruleList().size())
	  {
	    // At this point we have the fragment oligomer. However, do
	    // not forget that the user might ask for fragments that
	    // bear more than the single charge that was intrinsically
	    // computed within the formula of the fragmentation
	    // specification. 
      
	    // So, first create an oligomer with the "default"
	    // fragmentation specification-driven 1+ charge.  

	    FragmentOligomer *newOligomer = 
	      new FragmentOligomer(*oligomer1);
	  
	    QString name = QString("%1#%2#z=1")
	      .arg(fragOptions.name())
	      .arg(number + 1);
	  
	    newOligomer->setName(name);
	  
	    oligomerList->append(newOligomer);
	  
	    ++count;
	  
	    // And then create as many new oligomers as needed for the
	    // different charge levels asked by the user.  Because the
	    // ionization changes the values in the oligomer, and we
	    // need a new oligomer each time, we duplicate the oligomer
	    // each time we need it.

	    int startIonizeLevel = fragOptions.startIonizeLevel();
	    int endIonizeLevel = fragOptions.endIonizeLevel();

// 		qDebug() << __FILE__ << __LINE__ 
// 			 << "ionizeRule formula:" 
// 			 << ionizeRule.formula();

	    for (int kter = startIonizeLevel; kter < endIonizeLevel; ++kter)
	      {
		IonizeRule ionizeRule(m_ionizeRule);
		ionizeRule.setLevel(kter);

		FragmentOligomer *newOligomer = 
		  new FragmentOligomer(*oligomer1);
	  
		// If the result of the call below is -1, then that
		// means that there was an error and we should return
		// immediately. If it is 0, then that means that no
		// error was encountered, but that no actual ionization
		// took place, so we need not take into account the
		// oligomer.

		int res = newOligomer->ionize(ionizeRule);
	      
		if(res == -1)
		  {
		    delete newOligomer;
		    delete oligomer1;
		  
		    return -1;
		  }
		else if (res == 0)
		  {
		    delete newOligomer;
		
		    continue;
		  }
	      
		// At this point the ionization did indeed perform
		// something interesting, continue constructing the
		// resulting oligomer.

		QString name = QString("%1#%2#z=%3")
		  .arg(fragOptions.name())
		  .arg(number + 1)
		  .arg(newOligomer->charge());
	  
		newOligomer->setName(name);
	  
		oligomerList->append(newOligomer);
	  
		++count;
	      }

	    // We can delete the template oligomer.
	    delete oligomer1;
	  }
	else
	  {
	    if (success)
	      delete oligomer1;
	    else
	      {
		// At this point we have the fragment oligomer. However, do
		// not forget that the user might ask for fragments that
		// bear more than the single charge that was intrinsically
		// computed within the formula of the fragmentation
		// specification. 

		// So, first create an oligomer with the "default"
		// fragmentation specification-driven 1+ charge.  

		FragmentOligomer *newOligomer = 
		  new FragmentOligomer(*oligomer1);
	  
		QString name = QString("%1%2#z=1")
		  .arg(fragOptions.name())
		  .arg(number + 1);
	  
		newOligomer->setName(name);
	  
		oligomerList->append(newOligomer);
	  
		++count;
	  
		// And then create as many new oligomers as needed for the
		// different charge levels asked by the user.  Because the
		// ionization changes the values in the oligomer, and we
		// need a new oligomer each time, we duplicate the oligomer
		// each time we need it.

		int startIonizeLevel = fragOptions.startIonizeLevel();
		int endIonizeLevel = fragOptions.endIonizeLevel();

		for(int kter = startIonizeLevel; kter < endIonizeLevel; ++kter)
		  {
		    IonizeRule ionizeRule(m_ionizeRule);
		    ionizeRule.setLevel(kter);

		    FragmentOligomer *newOligomer = 
		      new FragmentOligomer(*oligomer1);
	  
		    // If the result of the call below is -1, then that
		    // means that there was an error and we should return
		    // immediately. If it is 0, then that means that no
		    // error was encountered, but that no actual ionization
		    // took place, so we need not take into account the
		    // oligomer.

		    int res = newOligomer->ionize(ionizeRule);
	      
		    if (res == -1)
		      {
			delete newOligomer;
			delete oligomer1;
		  
			return -1;
		      }
		    else if (res == 0)
		      {
			delete newOligomer;
		
			continue;
		      }
	      
		    // At this point the ionization did indeed perform
		    // something interesting, continue constructing the
		    // resulting oligomer.

		    QString name = QString("%1%2#z=%3")
		      .arg(fragOptions.name())
		      .arg(number + 1)
		      .arg(newOligomer->charge());
	  
		    newOligomer->setName(name);
	  
		    oligomerList->append(newOligomer);
	  
		    ++count;
		  }

		// We can delete the template oligomer.
		delete oligomer1;
	      }
	  }
      }
    // End of
    //  for (int iter = fragOptions.endIndex(); 
    //  iter > fragOptions.endIndex() - 1; --iter, ++number)
  
    m_oligomerLists.append(oligomerList);

    return count;
  }


  bool 
  Fragmenter::accountFragRule(FragRule *fragRule, bool onlyCheck,
			       int index, int fragEnd,
			       Ponderable *ponderable)
  {
    const Monomer *prevMonomer = 0;
    const Monomer *nextMonomer = 0;
  
    QList<Atom *> refList = mp_polChemDef->atomList();

    Q_ASSERT(fragRule);
  
    if (!onlyCheck)
      Q_ASSERT(ponderable);
    

    const Monomer *monomer = mp_polymer->at(index);
  
    if (!fragRule->currCode().isEmpty())
      if (fragRule->currCode() != monomer->code())
	return false;
  
    if (!fragRule->prevCode().isEmpty() && 
	!fragRule->nextCode().isEmpty())
      {
	if(fragEnd & MXT_FRAG_END_LEFT || fragEnd & MXT_FRAG_END_NONE)
	  {
	    if (!index)
	      // There cannot be any prevCode since we are at index ==
	      // 0, at the first monomer of the fragmentation
	      // series. That means that we can return immediately.
	      return false;

	    // Since we know that we are either in LEFT or NONE end
	    // mode, we know that previous is at index 'index' - 1. Thus
	    // get the monomer out of the sequence for this index.

	    prevMonomer = mp_polymer->at(index - 1);
	  
	    if (index == mp_polymer->size() - 1)
	      // There cannot be any next code since we are already at
	      // the last monomer in the fragmentation series.
	      return false;

	    nextMonomer = mp_polymer->at(index + 1);
	  }
	else if (fragEnd & MXT_FRAG_END_RIGHT)
	  {
	    if (!index)
	      // There cannot be any nextCode since currCode is the last
	      // monomer in the fragmentation series.
	      return false;
	  
	    nextMonomer = mp_polymer->at(index - 1);
	  
	    if (index == mp_polymer->size() - 1)
	      // There cannot be any previous code since currCode is the
	      // first in the fragmentation series.
	      return false;
	  
	    prevMonomer = mp_polymer->at(index + 1);
	  }
	else
	  return false;
      
	// Now that the prevCode and nextCode have been correctly
	// identified, we can go on and check if some conditions are
	// met.

	if(fragRule->prevCode() == prevMonomer->code() &&
	    fragRule->nextCode() == nextMonomer->code())
	  {
	    if (onlyCheck)
	      return true;
	  
	    // The fragmentation rule condition is met, we can apply its
	    // formula.

	    if (!fragRule->Formula::accountMasses(refList,
						   ponderable))
	      {
		qDebug() << __FILE__ << __LINE__
			  << "Failed to account fragmentation rule";
	      
		return false;
	      }
	  
	    return true;
	  }
	else
	  {
	    if (onlyCheck)
	      return false;
	    else
	      return true;
	  }
      }
    // End of 
    //   if (!fragRule->prevCode().isEmpty() && 
    //   !fragRule->nextCode().isEmpty())
    else if (!fragRule->prevCode().isEmpty())
      {
	if(fragEnd & MXT_FRAG_END_LEFT || fragEnd & MXT_FRAG_END_NONE)
	  {
	    if (!index)
	      // There cannot be any prevCode since currCode is already
	      // the first of the fragmentation series.
	      return false;
	  
	    // Since we know that fragEnd is either LEFT or NONE end, we
	    // know what index has the prevCode:
	  
	    prevMonomer = mp_polymer->at(index - 1);
	  }
	else if (fragEnd & MXT_FRAG_END_RIGHT)
	  {
	    if (index == mp_polymer->size() - 1)
	      // There cannot be any prevCode since currCode is already
	      // the first of the fragmentation series.
	      return false;
	    
	    prevMonomer = mp_polymer->at(index + 1);
	  }
	else
	  return false;
  
	// Now that we have correctly identified the prevCode, we can go
	// on and check if some conditions are met.
      
	if(fragRule->prevCode() == prevMonomer->code())
	  {
	    if (onlyCheck)
	      return true;
	  
	    // The fragmentation rule condition is met, we can apply its
	    // formula.

	    if (!fragRule->Formula::accountMasses(refList, ponderable))
	      {
		qDebug() << __FILE__ << __LINE__
			  << "Failed to account fragmentation rule";
	      
		return false;
	      } 
	  
	    return true;
	  }
	else
	  {
	    if (onlyCheck)
	      return false;
	    else
	      return true;
	  }
      }
    // End of 
    // else if (!fragRule->prevCode().isEmpty())
    else if (!fragRule->nextCode().isEmpty())
      {
	if(fragEnd & MXT_FRAG_END_LEFT || fragEnd & MXT_FRAG_END_NONE)
	  {
	    if (index == mp_polymer->size() - 1)
	      // There cannot be any nextCode since currCode is already
	      // the last of the fragmentation series.
	      return false;
	  
	    // Since we know that fragEnd is either LEFT or NONE end, we
	    // know what index has the prevCode:
	  
	    nextMonomer = mp_polymer->at(index + 1);
	  }
	else if (fragEnd & MXT_FRAG_END_RIGHT)
	  {
	    if (!index)
	      // There cannot be any prevCode since currCode is already
	      // the last of the fragmentation series.
	      return false;
	    
	    nextMonomer = mp_polymer->at(index - 1);
	  }
	else
	  return false;
  
	// Now that we have correctly identified the nextCode, we can go
	// on and check if some conditions are met.
      
	if(fragRule->nextCode() == nextMonomer->code())
	  {
	    if (onlyCheck)
	      return true;
	  
	    // The fragmentation rule condition is met, we can apply its
	    // formula.

	    if (!fragRule->Formula::accountMasses(refList,
						   ponderable))
	      {
		qDebug() << __FILE__ << __LINE__
			  << "Failed to account fragmentation rule";
	      
		return false;
	      }
	  
	    return true;
	  }
	else
	  {
	    if (onlyCheck)
	      return false;
	    else
	      return true;
	  }
      }
    // End of 
    // else if (!fragRule->nextCode().isEmpty())
    else
      {
	// All the prev and next codes are empty, which means that we
	// consider the conditions verified.
	if(onlyCheck)
	  return true;
      
	if(!fragRule->Formula::accountMasses(refList,
					       ponderable))
	  {
	    qDebug() << __FILE__ << __LINE__
		      << "Failed to account fragmentation rule";
	  
	    return false;
	  }
      
	return true;
      }
      
    // We should never reach this point !
    Q_ASSERT(0);
  
    return false;
  }




  void 
  Fragmenter::freeOligomerLists()
  {
    // Iterate in the QList of OligomerList objects and free each one
    // in turn.

    while(!m_oligomerLists.isEmpty())
      {
	OligomerList *currentList = m_oligomerLists.takeFirst();
	Q_ASSERT(currentList);
      
	delete currentList;
      }
  }

} // namespace massXpert