// @(#)root/matrix:$Id$
// Authors: Fons Rademakers, Eddy Offermann  Nov 2003

/*************************************************************************
 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
 * All rights reserved.                                                  *
 *                                                                       *
 * For the licensing terms see $ROOTSYS/LICENSE.                         *
 * For the list of contributors see $ROOTSYS/README/CREDITS.             *
 *************************************************************************/

//////////////////////////////////////////////////////////////////////////
//                                                                      //
// TVectorT                                                             //
//                                                                      //
// Template class of Vectors in the linear algebra package              //
//                                                                      //
// Unless otherwise specified, vector indices always start with 0,      //
// spanning up to the specified limit-1.                                //
//                                                                      //
// For (n) vectors where n <= kSizeMax (5 currently) storage space is   //
// available on the stack, thus avoiding expensive allocation/          //
// deallocation of heap space . However, this introduces of course      //
// kSizeMax overhead for each vector object . If this is an issue       //
// recompile with a new appropriate value (>=0) for kSizeMax            //
//                                                                      //
// Another way to assign and store vector data is through Use           //
// see for instance stressLinear.cxx file .                             //
//                                                                      //
// Note that Constructors/assignments exists for all different matrix   //
// views                                                                //
//                                                                      //
// For usage examples see $ROOTSYS/test/stressLinear.cxx                //
//                                                                      //
//////////////////////////////////////////////////////////////////////////

#include "TVectorT.h"
#include "TClass.h"
#include "TMath.h"
#include "TROOT.h"
#include "Varargs.h"

templateClassImp(TVectorT)


//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::Delete_m(Int_t size,Element *&m)
{
// Delete data pointer m, if it was assigned on the heap

   if (m) {
      if (size > kSizeMax)
         delete [] m;
      m = 0;
   }
}

//______________________________________________________________________________
template<class Element>
Element* TVectorT<Element>::New_m(Int_t size)
{
// Return data pointer . if requested size <= kSizeMax, assign pointer
// to the stack space

   if (size == 0) return 0;
   else {
      if ( size <= kSizeMax )
         return fDataStack;
      else {
         Element *heap = new Element[size];
         return heap;
      }
   }
}

//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::Add(const TVectorT<Element> &v)
{
// Add vector v to this vector

   if (gMatrixCheck && !AreCompatible(*this,v)) {
      Error("Add(TVectorT<Element> &)","vector's not compatible");
      return;
   }

   const Element *sp = v.GetMatrixArray();
         Element *tp = this->GetMatrixArray();
   const Element * const tp_last = tp+fNrows;
   while (tp < tp_last)
      *tp++ += *sp++;
}

//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::Add(const TVectorT<Element> &v1,const TVectorT<Element> &v2)
{
// Set this vector to v1+v2

   if (gMatrixCheck) {
      if ( !AreCompatible(*this,v1) && !AreCompatible(*this,v2)) {
         Error("Add(TVectorT<Element> &)","vectors not compatible");
         return;
      }
   }

   const Element *sv1 = v1.GetMatrixArray();
   const Element *sv2 = v2.GetMatrixArray();
         Element *tp = this->GetMatrixArray();
   const Element * const tp_last = tp+fNrows;
   while (tp < tp_last)
      *tp++ = *sv1++ + *sv2++;
}

//______________________________________________________________________________
template<class Element>
Int_t TVectorT<Element>::Memcpy_m(Element *newp,const Element *oldp,Int_t copySize,
                                  Int_t newSize,Int_t oldSize)
{
// Copy copySize doubles from *oldp to *newp . However take care of the
// situation where both pointers are assigned to the same stack space

   if (copySize == 0 || oldp == newp) return 0;
   else {
      if ( newSize <= kSizeMax && oldSize <= kSizeMax ) {
         // both pointers are inside fDataStack, be careful with copy direction !
         if (newp > oldp) {
            for (Int_t i = copySize-1; i >= 0; i--)
               newp[i] = oldp[i];
         } else {
            for (Int_t i = 0; i < copySize; i++)
               newp[i] = oldp[i];
         }
      }
      else {
         memcpy(newp,oldp,copySize*sizeof(Element));
      }
   }
   return 0;
}

//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::Allocate(Int_t nrows,Int_t row_lwb,Int_t init)
{
// Allocate new vector. Arguments are number of rows and row
// lowerbound (0 default).

   fIsOwner  = kTRUE;
   fNrows    = 0;
   fRowLwb   = 0;
   fElements = 0;

   if (nrows < 0)
   {
      Error("Allocate","nrows=%d",nrows);
      return;
   }

   MakeValid();
   fNrows   = nrows;
   fRowLwb  = row_lwb;

   fElements = New_m(fNrows);
   if (init)
      memset(fElements,0,fNrows*sizeof(Element));
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(Int_t n)
{
// Constructor n-vector

   Allocate(n,0,1);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(Int_t lwb,Int_t upb)
{
// Constructor [lwb..upb]-vector

   Allocate(upb-lwb+1,lwb,1);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(Int_t n,const Element *elements)
{
// Constructor n-vector with data copied from array elements

   Allocate(n,0);
   SetElements(elements);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(Int_t lwb,Int_t upb,const Element *elements)
{
// Constructor [lwb..upb]-vector with data copied from array elements

   Allocate(upb-lwb+1,lwb);
   SetElements(elements);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(const TVectorT &another) : TObject(another)
{
// Copy constructor

   R__ASSERT(another.IsValid());
   Allocate(another.GetUpb()-another.GetLwb()+1,another.GetLwb());
   *this = another;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(const TMatrixTRow_const<Element> &mr) : TObject()
{
// Constructor : create vector from matrix row

   const TMatrixTBase<Element> *mt = mr.GetMatrix();
   R__ASSERT(mt->IsValid());
   Allocate(mt->GetColUpb()-mt->GetColLwb()+1,mt->GetColLwb());
   *this = mr;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(const TMatrixTColumn_const<Element> &mc) : TObject()
{
// Constructor : create vector from matrix column

   const TMatrixTBase<Element> *mt = mc.GetMatrix();
   R__ASSERT(mt->IsValid());
   Allocate(mt->GetRowUpb()-mt->GetRowLwb()+1,mt->GetRowLwb());
   *this = mc;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(const TMatrixTDiag_const<Element> &md) : TObject()
{
// Constructor : create vector from matrix diagonal

   const TMatrixTBase<Element> *mt = md.GetMatrix();
   R__ASSERT(mt->IsValid());
   Allocate(TMath::Min(mt->GetNrows(),mt->GetNcols()));
   *this = md;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element>::TVectorT(Int_t lwb,Int_t upb,Element va_(iv1), ...)
{
// Make a vector and assign initial values. Argument list should contain
// Element values to assign to vector elements. The list must be
// terminated by the string "END". Example:
// TVectorT foo(1,3,0.0,1.0,1.5,"END");

   if (upb < lwb) {
      Error("TVectorT(Int_t, Int_t, ...)","upb(%d) < lwb(%d)",upb,lwb);
      return;
   }

   Allocate(upb-lwb+1,lwb);

   va_list args;
   va_start(args,va_(iv1));             // Init 'args' to the beginning of
                                        // the variable length list of args

   (*this)(lwb) = iv1;
   for (Int_t i = lwb+1; i <= upb; i++)
      (*this)(i) = (Element)va_arg(args,Double_t);

   if (strcmp((char *)va_arg(args,char *),"END"))
      Error("TVectorT(Int_t, Int_t, ...)","argument list must be terminated by \"END\"");

   va_end(args);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::ResizeTo(Int_t lwb,Int_t upb)
{
// Resize the vector to [lwb:upb] .
// New dynamic elemenst are created, the overlapping part of the old ones are
// copied to the new structures, then the old elements are deleted.

   R__ASSERT(IsValid());
   if (!fIsOwner) {
      Error("ResizeTo(lwb,upb)","Not owner of data array,cannot resize");
      return *this;
   }

   const Int_t new_nrows = upb-lwb+1;

   if (fNrows > 0) {

      if (fNrows == new_nrows && fRowLwb == lwb)
         return *this;
      else if (new_nrows == 0) {
         Clear();
         return *this;
      }

      Element    *elements_old = GetMatrixArray();
      const Int_t  nrows_old    = fNrows;
      const Int_t  rowLwb_old   = fRowLwb;

      Allocate(new_nrows,lwb);
      R__ASSERT(IsValid());
      if (fNrows > kSizeMax || nrows_old > kSizeMax)
         memset(GetMatrixArray(),0,fNrows*sizeof(Element));
      else if (fNrows > nrows_old)
         memset(GetMatrixArray()+nrows_old,0,(fNrows-nrows_old)*sizeof(Element));

    // Copy overlap
      const Int_t rowLwb_copy = TMath::Max(fRowLwb,rowLwb_old);
      const Int_t rowUpb_copy = TMath::Min(fRowLwb+fNrows-1,rowLwb_old+nrows_old-1);
      const Int_t nrows_copy  = rowUpb_copy-rowLwb_copy+1;

      const Int_t nelems_new = fNrows;
      Element *elements_new = GetMatrixArray();
      if (nrows_copy > 0) {
         const Int_t rowOldOff = rowLwb_copy-rowLwb_old;
         const Int_t rowNewOff = rowLwb_copy-fRowLwb;
         Memcpy_m(elements_new+rowNewOff,elements_old+rowOldOff,nrows_copy,nelems_new,nrows_old);
      }

      Delete_m(nrows_old,elements_old);
   } else {
      Allocate(upb-lwb+1,lwb,1);
   }

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::Use(Int_t lwb,Int_t upb,Element *data)
{
// Use the array data to fill the vector lwb..upb]

   if (upb < lwb) {
      Error("Use","upb(%d) < lwb(%d)",upb,lwb);
      return *this;
   }

   Clear();
   fNrows    = upb-lwb+1;
   fRowLwb   = lwb;
   fElements = data;
   fIsOwner  = kFALSE;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::GetSub(Int_t row_lwb,Int_t row_upb,TVectorT<Element> &target,Option_t *option) const
{
// Get subvector [row_lwb..row_upb]; The indexing range of the
// returned vector depends on the argument option:
//
// option == "S" : return [0..row_upb-row_lwb+1] (default)
// else          : return [row_lwb..row_upb]

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      if (row_lwb < fRowLwb || row_lwb > fRowLwb+fNrows-1) {
         Error("GetSub","row_lwb out of bounds");
         return target;
      }
      if (row_upb < fRowLwb || row_upb > fRowLwb+fNrows-1) {
         Error("GetSub","row_upb out of bounds");
         return target;
      }
      if (row_upb < row_lwb) {
         Error("GetSub","row_upb < row_lwb");
         return target;
      }
   }

   TString opt(option);
   opt.ToUpper();
   const Int_t shift = (opt.Contains("S")) ? 1 : 0;

   Int_t row_lwb_sub;
   Int_t row_upb_sub;
   if (shift) {
      row_lwb_sub = 0;
      row_upb_sub = row_upb-row_lwb;
   } else {
      row_lwb_sub = row_lwb;
      row_upb_sub = row_upb;
   }

   target.ResizeTo(row_lwb_sub,row_upb_sub);
   const Int_t nrows_sub = row_upb_sub-row_lwb_sub+1;

   const Element *ap = this->GetMatrixArray()+(row_lwb-fRowLwb);
         Element *bp = target.GetMatrixArray();

   for (Int_t irow = 0; irow < nrows_sub; irow++)
      *bp++ = *ap++;

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::SetSub(Int_t row_lwb,const TVectorT<Element> &source)
{
// Insert vector source starting at [row_lwb], thereby overwriting the part
// [row_lwb..row_lwb+nrows_source];

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      R__ASSERT(source.IsValid());

      if (row_lwb < fRowLwb && row_lwb > fRowLwb+fNrows-1) {
         Error("SetSub","row_lwb outof bounds");
         return *this;
      }
      if (row_lwb+source.GetNrows() > fRowLwb+fNrows) {
         Error("SetSub","source vector too large");
         return *this;
      }
   }

   const Int_t nRows_source = source.GetNrows();

   const Element *bp = source.GetMatrixArray();
         Element *ap = this->GetMatrixArray()+(row_lwb-fRowLwb);

   for (Int_t irow = 0; irow < nRows_source; irow++)
      *ap++ = *bp++;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::Zero()
{
// Set vector elements to zero

   R__ASSERT(IsValid());
   memset(this->GetMatrixArray(),0,fNrows*sizeof(Element));
   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::Abs()
{
// Take an absolute value of a vector, i.e. apply Abs() to each element.

   R__ASSERT(IsValid());

         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      *ep = TMath::Abs(*ep);
      ep++;
   }

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::Sqr()
{
// Square each element of the vector.

   R__ASSERT(IsValid());

         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      *ep = (*ep) * (*ep);
      ep++;
   }

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::Sqrt()
{
// Take square root of all elements.

   R__ASSERT(IsValid());

         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      R__ASSERT(*ep >= 0);
      if (*ep >= 0)
         *ep = TMath::Sqrt(*ep);
      else
         Error("Sqrt()","v(%ld) = %g < 0",Long_t(ep-this->GetMatrixArray()),(float)*ep);
      ep++;
   }

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::Invert()
{
// v[i] = 1/v[i]

   R__ASSERT(IsValid());

         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      R__ASSERT(*ep != 0.0);
      if (*ep != 0.0)
         *ep = 1./ *ep;
      else
         Error("Invert()","v(%ld) = %g",Long_t(ep-this->GetMatrixArray()),(float)*ep);
      ep++;
   }

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::SelectNonZeros(const TVectorT<Element> &select)
{
// Keep only element as selected through array select non-zero

   if (gMatrixCheck && !AreCompatible(*this,select)) {
      Error("SelectNonZeros(const TVectorT<Element> &","vector's not compatible");
      return *this;
   }

   const Element *sp = select.GetMatrixArray();
         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      if (*sp == 0.0)
         *ep = 0.0;
      sp++; ep++;
   }

   return *this;
}

//______________________________________________________________________________
template<class Element>
Element TVectorT<Element>::Norm1() const
{
// Compute the 1-norm of the vector SUM{ |v[i]| }.

   R__ASSERT(IsValid());

   Element norm = 0;
   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      norm += TMath::Abs(*ep++);

   return norm;
}

//______________________________________________________________________________
template<class Element>
Element TVectorT<Element>::Norm2Sqr() const
{
// Compute the square of the 2-norm SUM{ v[i]^2 }.

   R__ASSERT(IsValid());

   Element norm = 0;
   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      norm += (*ep) * (*ep);
      ep++;
   }

   return norm;
}

//______________________________________________________________________________
template<class Element>
Element TVectorT<Element>::NormInf() const
{
// Compute the infinity-norm of the vector MAX{ |v[i]| }.

   R__ASSERT(IsValid());

   Element norm = 0;
   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      norm = TMath::Max(norm,TMath::Abs(*ep++));

   return norm;
}

//______________________________________________________________________________
template<class Element>
Int_t TVectorT<Element>::NonZeros() const
{
// Compute the number of elements != 0.0

   R__ASSERT(IsValid());

   Int_t nr_nonzeros = 0;
   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      if (*ep++) nr_nonzeros++;

   return nr_nonzeros;
}

//______________________________________________________________________________
template<class Element>
Element TVectorT<Element>::Sum() const
{
// Compute sum of elements

   R__ASSERT(IsValid());

   Element sum = 0.0;
   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      sum += *ep++;

   return sum;
}

//______________________________________________________________________________
template<class Element>
Element TVectorT<Element>::Min() const
{
// return minimum vector element value

   R__ASSERT(IsValid());

   const Int_t index = TMath::LocMin(fNrows,fElements);
   return fElements[index];
}

//______________________________________________________________________________
template<class Element>
Element TVectorT<Element>::Max() const
{
// return maximum vector element value

   R__ASSERT(IsValid());

   const Int_t index = TMath::LocMax(fNrows,fElements);
   return fElements[index];
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator=(const TVectorT<Element> &source)
{
// Notice that this assignment does NOT change the ownership :
// if the storage space was adopted, source is copied to
// this space .

   if (gMatrixCheck && !AreCompatible(*this,source)) {
      Error("operator=(const TVectorT<Element> &)","vectors not compatible");
      return *this;
   }

   if (this->GetMatrixArray() != source.GetMatrixArray()) {
      TObject::operator=(source);
      memcpy(fElements,source.GetMatrixArray(),fNrows*sizeof(Element));
   }
   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator=(const TMatrixTRow_const<Element> &mr)
{
// Assign a matrix row to a vector.

   const TMatrixTBase<Element> *mt = mr.GetMatrix();

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      R__ASSERT(mt->IsValid());
      if (mt->GetColLwb() != fRowLwb || mt->GetNcols() != fNrows) {
         Error("operator=(const TMatrixTRow_const &)","vector and row not compatible");
         return *this;
      }
   }

   const Int_t inc   = mr.GetInc();
   const Element *rp = mr.GetPtr();              // Row ptr
         Element *ep = this->GetMatrixArray();   // Vector ptr
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      *ep++ = *rp;
       rp += inc;
   }

   R__ASSERT(rp == mr.GetPtr()+mt->GetNcols());

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator=(const TMatrixTColumn_const<Element> &mc)
{
// Assign a matrix column to a vector.

   const TMatrixTBase<Element> *mt = mc.GetMatrix();

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      R__ASSERT(mt->IsValid());
      if (mt->GetRowLwb() != fRowLwb || mt->GetNrows() != fNrows) {
         Error("operator=(const TMatrixTColumn_const &)","vector and column not compatible");
         return *this;
      }
   }

   const Int_t inc    = mc.GetInc();
   const Element *cp = mc.GetPtr();              // Column ptr
         Element *ep = this->GetMatrixArray();   // Vector ptr
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      *ep++ = *cp;
       cp += inc;
   }

   R__ASSERT(cp == mc.GetPtr()+mt->GetNoElements());

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator=(const TMatrixTDiag_const<Element> &md)
{
// Assign the matrix diagonal to a vector.

   const TMatrixTBase<Element> *mt = md.GetMatrix();

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      R__ASSERT(mt->IsValid());
      if (md.GetNdiags() != fNrows) {
         Error("operator=(const TMatrixTDiag_const &)","vector and matrix-diagonal not compatible");
        return *this;
      }
   }

   const Int_t    inc = md.GetInc();
   const Element *dp = md.GetPtr();              // Diag ptr
         Element *ep = this->GetMatrixArray();   // Vector ptr
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      *ep++ = *dp;
       dp += inc;
   }

   R__ASSERT(dp < md.GetPtr()+mt->GetNoElements()+inc);

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator=(const TMatrixTSparseRow_const<Element> &mr)
{
// Assign a sparse matrix row to a vector. The matrix row is implicitly transposed
// to allow the assignment in the strict sense.

   const TMatrixTBase<Element> *mt = mr.GetMatrix();

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      R__ASSERT(mt->IsValid());
      if (mt->GetColLwb() != fRowLwb || mt->GetNcols() != fNrows) {
         Error("operator=(const TMatrixTSparseRow_const &)","vector and row not compatible");
         return *this;
      }
   }

   const Int_t nIndex = mr.GetNindex();
   const Element * const prData = mr.GetDataPtr();          // Row Data ptr
   const Int_t    * const prCol  = mr.GetColPtr();           // Col ptr
         Element * const pvData = this->GetMatrixArray();   // Vector ptr

   memset(pvData,0,fNrows*sizeof(Element));
   for (Int_t index = 0; index < nIndex; index++) {
      const Int_t icol = prCol[index];
      pvData[icol] = prData[index];
   }

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator=(const TMatrixTSparseDiag_const<Element> &md)
{
// Assign a sparse matrix diagonal to a vector.

  const TMatrixTBase<Element> *mt = md.GetMatrix();

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      R__ASSERT(mt->IsValid());
      if (md.GetNdiags() != fNrows) {
         Error("operator=(const TMatrixTSparseDiag_const &)","vector and matrix-diagonal not compatible");
         return *this;
      }
   }

   Element * const pvData = this->GetMatrixArray();
   for (Int_t idiag = 0; idiag < fNrows; idiag++)
      pvData[idiag] = md(idiag);

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator=(Element val)
{
// Assign val to every element of the vector.

   R__ASSERT(IsValid());

         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      *ep++ = val;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator+=(Element val)
{
// Add val to every element of the vector.

   R__ASSERT(IsValid());

         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      *ep++ += val;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator-=(Element val)
{
// Subtract val from every element of the vector.

   R__ASSERT(IsValid());

         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      *ep++ -= val;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator*=(Element val)
{
// Multiply every element of the vector with val.

   R__ASSERT(IsValid());

         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      *ep++ *= val;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator+=(const TVectorT<Element> &source)
{
// Add vector source

   if (gMatrixCheck && !AreCompatible(*this,source)) {
      Error("operator+=(const TVectorT<Element> &)","vector's not compatible");
      return *this;
   }

   const Element *sp = source.GetMatrixArray();
         Element *tp = this->GetMatrixArray();
   const Element * const tp_last = tp+fNrows;
   while (tp < tp_last)
      *tp++ += *sp++;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator-=(const TVectorT<Element> &source)
{
// Subtract vector source

   if (gMatrixCheck && !AreCompatible(*this,source)) {
      Error("operator-=(const TVectorT<Element> &)","vector's not compatible");
      return *this;
   }

   const Element *sp = source.GetMatrixArray();
         Element *tp = this->GetMatrixArray();
   const Element * const tp_last = tp+fNrows;
   while (tp < tp_last)
      *tp++ -= *sp++;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator*=(const TMatrixT<Element> &a)
{
// "Inplace" multiplication target = A*target. A needn't be a square one
// If target has to be resized, it should own the storage: fIsOwner = kTRUE

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      R__ASSERT(a.IsValid());
      if (a.GetNcols() != fNrows || a.GetColLwb() != fRowLwb) {
         Error("operator*=(const TMatrixT &)","vector and matrix incompatible");
         return *this;
      }
   }

   const Bool_t doResize = (fNrows != a.GetNrows() || fRowLwb != a.GetRowLwb());
   if (doResize && !fIsOwner) {
      Error("operator*=(const TMatrixT &)","vector has to be resized but not owner");
      return *this;
   }

   Element work[kWorkMax];
   Bool_t isAllocated = kFALSE;
   Element *elements_old = work;
   const Int_t nrows_old = fNrows;
   if (nrows_old > kWorkMax) {
      isAllocated = kTRUE;
      elements_old = new Element[nrows_old];
   }
   memcpy(elements_old,fElements,nrows_old*sizeof(Element));

   if (doResize) {
      const Int_t rowlwb_new = a.GetRowLwb();
      const Int_t nrows_new  = a.GetNrows();
      ResizeTo(rowlwb_new,rowlwb_new+nrows_new-1);
   }
   memset(fElements,0,fNrows*sizeof(Element));

   const Element *mp = a.GetMatrixArray();     // Matrix row ptr
         Element *tp = this->GetMatrixArray(); // Target vector ptr
#ifdef CBLAS
   if (typeid(Element) == typeid(Double_t))
      cblas_dgemv(CblasRowMajor,CblasNoTrans,a.GetNrows(),a.GetNcols(),1.0,mp,
                  a.GetNcols(),elements_old,1,0.0,tp,1);
   else if (typeid(Element) != typeid(Float_t))
      cblas_sgemv(CblasRowMajor,CblasNoTrans,a.GetNrows(),a.GetNcols(),1.0,mp,
                  a.GetNcols(),elements_old,1,0.0,tp,1);
   else
      Error("operator*=","type %s not implemented in BLAS library",typeid(Element));
#else
   const Element * const tp_last = tp+fNrows;
   while (tp < tp_last) {
      Element sum = 0;
      for (const Element *sp = elements_old; sp < elements_old+nrows_old; )
         sum += *sp++ * *mp++;
      *tp++ = sum;
   }
   R__ASSERT(mp == a.GetMatrixArray()+a.GetNoElements());
#endif

   if (isAllocated)
      delete [] elements_old;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator*=(const TMatrixTSparse<Element> &a)
{
// "Inplace" multiplication target = A*target. A needn't be a square one
// If target has to be resized, it should own the storage: fIsOwner = kTRUE

   if (gMatrixCheck) {
      R__ASSERT(IsValid());
      R__ASSERT(a.IsValid());
      if (a.GetNcols() != fNrows || a.GetColLwb() != fRowLwb) {
         Error("operator*=(const TMatrixTSparse &)","vector and matrix incompatible");
         return *this;
      }
   }

   const Bool_t doResize = (fNrows != a.GetNrows() || fRowLwb != a.GetRowLwb());
   if (doResize && !fIsOwner) {
      Error("operator*=(const TMatrixTSparse &)","vector has to be resized but not owner");
      return *this;
   }

   Element work[kWorkMax];
   Bool_t isAllocated = kFALSE;
   Element *elements_old = work;
   const Int_t nrows_old = fNrows;
   if (nrows_old > kWorkMax) {
      isAllocated = kTRUE;
      elements_old = new Element[nrows_old];
   }
   memcpy(elements_old,fElements,nrows_old*sizeof(Element));

   if (doResize) {
      const Int_t rowlwb_new = a.GetRowLwb();
      const Int_t nrows_new  = a.GetNrows();
      ResizeTo(rowlwb_new,rowlwb_new+nrows_new-1);
   }
   memset(fElements,0,fNrows*sizeof(Element));

   const Int_t   * const pRowIndex = a.GetRowIndexArray();
   const Int_t   * const pColIndex = a.GetColIndexArray();
   const Element * const mp        = a.GetMatrixArray();     // Matrix row ptr

   const Element * const sp = elements_old;
         Element *       tp = this->GetMatrixArray(); // Target vector ptr

   for (Int_t irow = 0; irow < fNrows; irow++) {
      const Int_t sIndex = pRowIndex[irow];
      const Int_t eIndex = pRowIndex[irow+1];
      Element sum = 0.0;
      for (Int_t index = sIndex; index < eIndex; index++) {
         const Int_t icol = pColIndex[index];
         sum += mp[index]*sp[icol];
      }
      tp[irow] = sum;
   }

   if (isAllocated)
      delete [] elements_old;

   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::operator*=(const TMatrixTSym<Element> &a)
{
// "Inplace" multiplication target = A*target. A is symmetric .
// vector size will not change

  if (gMatrixCheck) {
     R__ASSERT(IsValid());
     R__ASSERT(a.IsValid());
     if (a.GetNcols() != fNrows || a.GetColLwb() != fRowLwb) {
        Error("operator*=(const TMatrixTSym &)","vector and matrix incompatible");
        return *this;
     }
  }

   Element work[kWorkMax];
   Bool_t isAllocated = kFALSE;
   Element *elements_old = work;
   const Int_t nrows_old = fNrows;
   if (nrows_old > kWorkMax) {
      isAllocated = kTRUE;
      elements_old = new Element[nrows_old];
   }
   memcpy(elements_old,fElements,nrows_old*sizeof(Element));
   memset(fElements,0,fNrows*sizeof(Element));

   const Element *mp = a.GetMatrixArray();     // Matrix row ptr
         Element *tp = this->GetMatrixArray(); // Target vector ptr
#ifdef CBLAS
   if (typeid(Element) == typeid(Double_t))
      cblas_dsymv(CblasRowMajor,CblasUpper,fNrows,1.0,mp,
                  fNrows,elements_old,1,0.0,tp,1);
   else if (typeid(Element) != typeid(Float_t))
      cblas_ssymv(CblasRowMajor,CblasUpper,fNrows,1.0,mp,
                  fNrows,elements_old,1,0.0,tp,1);
   else
      Error("operator*=","type %s not implemented in BLAS library",typeid(Element));
#else
   const Element * const tp_last = tp+fNrows;
   while (tp < tp_last) {
      Element sum = 0;
      for (const Element *sp = elements_old; sp < elements_old+nrows_old; )
         sum += *sp++ * *mp++;
      *tp++ = sum;
   }
   R__ASSERT(mp == a.GetMatrixArray()+a.GetNoElements());
#endif

   if (isAllocated)
      delete [] elements_old;

   return *this;
}

//______________________________________________________________________________
template<class Element>
Bool_t TVectorT<Element>::operator==(Element val) const
{
// Are all vector elements equal to val?

   R__ASSERT(IsValid());

   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      if (!(*ep++ == val))
         return kFALSE;

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
Bool_t TVectorT<Element>::operator!=(Element val) const
{
// Are all vector elements not equal to val?

   R__ASSERT(IsValid());

   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      if (!(*ep++ != val))
         return kFALSE;

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
Bool_t TVectorT<Element>::operator<(Element val) const
{
// Are all vector elements < val?

   R__ASSERT(IsValid());

   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      if (!(*ep++ < val))
         return kFALSE;

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
Bool_t TVectorT<Element>::operator<=(Element val) const
{
// Are all vector elements <= val?

   R__ASSERT(IsValid());

   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      if (!(*ep++ <= val))
         return kFALSE;

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
Bool_t TVectorT<Element>::operator>(Element val) const
{
// Are all vector elements > val?

   R__ASSERT(IsValid());

   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      if (!(*ep++ > val))
         return kFALSE;

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
Bool_t TVectorT<Element>::operator>=(Element val) const
{
// Are all vector elements >= val?

   R__ASSERT(IsValid());

   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      if (!(*ep++ >= val))
         return kFALSE;

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
Bool_t TVectorT<Element>::MatchesNonZeroPattern(const TVectorT<Element> &select)
{
// Check if vector elements as selected through array select are non-zero

   if (gMatrixCheck && !AreCompatible(*this,select)) {
      Error("MatchesNonZeroPattern(const TVectorT&)","vector's not compatible");
      return kFALSE;
   }

   const Element *sp = select.GetMatrixArray();
   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      if (*sp == 0.0 && *ep != 0.0)
         return kFALSE;
      sp++; ep++;
   }

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
Bool_t TVectorT<Element>::SomePositive(const TVectorT<Element> &select)
{
// Check if vector elements as selected through array select are all positive

   if (gMatrixCheck && !AreCompatible(*this,select)) {
      Error("SomePositive(const TVectorT&)","vector's not compatible");
      return kFALSE;
   }

   const Element *sp = select.GetMatrixArray();
   const Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      if (*sp != 0.0 && *ep <= 0.0)
         return kFALSE;
      sp++; ep++;
   }

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::AddSomeConstant(Element val,const TVectorT<Element> &select)
{
// Add to vector elements as selected through array select the value val

   if (gMatrixCheck && !AreCompatible(*this,select))
      Error("AddSomeConstant(Element,const TVectorT&)(const TVectorT&)","vector's not compatible");

   const Element *sp = select.GetMatrixArray();
         Element *ep = this->GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp) {
      if (*sp)
         *ep += val;
      sp++; ep++;
   }
}

extern Double_t Drand(Double_t &ix);

//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::Randomize(Element alpha,Element beta,Double_t &seed)
{
// randomize vector elements value

   R__ASSERT(IsValid());

   const Element scale = beta-alpha;
   const Element shift = alpha/scale;

         Element *       ep = GetMatrixArray();
   const Element * const fp = ep+fNrows;
   while (ep < fp)
      *ep++ = scale*(Drand(seed)+shift);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::Apply(const TElementActionT<Element> &action)
{
// Apply action to each element of the vector.

   R__ASSERT(IsValid());
   for (Element *ep = fElements; ep < fElements+fNrows; ep++)
      action.Operation(*ep);
   return *this;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &TVectorT<Element>::Apply(const TElementPosActionT<Element> &action)
{
// Apply action to each element of the vector. In action the location
// of the current element is known.

   R__ASSERT(IsValid());

   Element *ep = fElements;
   for (action.fI = fRowLwb; action.fI < fRowLwb+fNrows; action.fI++)
      action.Operation(*ep++);

   R__ASSERT(ep == fElements+fNrows);

   return *this;
}

//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::Draw(Option_t *option)
{
// Draw this vector
// The histogram is named "TVectorT" by default and no title

   gROOT->ProcessLine(Form("THistPainter::PaintSpecialObjects((TObject*)0x%lx,\"%s\");",
                           (ULong_t)this, option));
}

//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::Print(Option_t *flag) const
{
// Print the vector as a list of elements.

  if (!IsValid()) {
      Error("Print","Vector is invalid");
      return;
   }

   printf("\nVector (%d) %s is as follows",fNrows,flag);

   printf("\n\n     |   %6d  |", 1);
   printf("\n%s\n", "------------------");
   for (Int_t i = 0; i < fNrows; i++) {
      printf("%4d |",i+fRowLwb);
      //printf("%11.4g \n",(*this)(i+fRowLwb));
      printf("%g \n",(*this)(i+fRowLwb));
   }
   printf("\n");
}

//______________________________________________________________________________
template<class Element>
Bool_t operator==(const TVectorT<Element> &v1,const TVectorT<Element> &v2)
{
// Check to see if two vectors are identical.

   if (!AreCompatible(v1,v2)) return kFALSE;
   return (memcmp(v1.GetMatrixArray(),v2.GetMatrixArray(),v1.GetNrows()*sizeof(Element)) == 0);
}

//______________________________________________________________________________
template<class Element>
Element operator*(const TVectorT<Element> &v1,const TVectorT<Element> &v2)
{
// Compute the scalar product.

   if (gMatrixCheck) {
      if (!AreCompatible(v1,v2)) {
         Error("operator*(const TVectorT<Element> &,const TVectorT<Element> &)","vector's are incompatible");
         return 0.0;
      }
   }

   return Dot(v1,v2);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> operator+(const TVectorT<Element> &source1,const TVectorT<Element> &source2)
{
// Return source1+source2

   TVectorT<Element> target = source1;
   target += source2;
   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> operator-(const TVectorT<Element> &source1,const TVectorT<Element> &source2)
{
// Return source1-source2

   TVectorT<Element> target = source1;
   target -= source2;
   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> operator*(const TMatrixT<Element> &a,const TVectorT<Element> &source)
{
// return A * source

   R__ASSERT(a.IsValid());
   TVectorT<Element> target(a.GetRowLwb(),a.GetRowUpb());
   return Add(target,Element(1.0),a,source);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> operator*(const TMatrixTSym<Element> &a,const TVectorT<Element> &source)
{
// return A * source

   R__ASSERT(a.IsValid());
   TVectorT<Element> target(a.GetRowLwb(),a.GetRowUpb());
   return Add(target,Element(1.0),a,source);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> operator*(const TMatrixTSparse<Element> &a,const TVectorT<Element> &source)
{
// return A * source

   R__ASSERT(a.IsValid());
   TVectorT<Element> target(a.GetRowLwb(),a.GetRowUpb());
   return Add(target,Element(1.0),a,source);
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> operator*(Element val,const TVectorT<Element> &source)
{
// return val * source

   TVectorT<Element> target = source;
   target *= val;
   return target;
}

//______________________________________________________________________________
template<class Element>
Element Dot(const TVectorT<Element> &v1,const TVectorT<Element> &v2)
{
// return inner-produvt v1 . v2

   const Element *v1p = v1.GetMatrixArray();
   const Element *v2p = v2.GetMatrixArray();
   Element sum = 0.0;
   const Element * const fv1p = v1p+v1.GetNrows();
   while (v1p < fv1p)
      sum += *v1p++ * *v2p++;

   return sum;
}

//______________________________________________________________________________
template <class Element1, class Element2>
TMatrixT<Element1>
OuterProduct(const TVectorT<Element1> &v1,const TVectorT<Element2> &v2)
{
   // Return the matrix M = v1 * v2'

   // TMatrixD::GetSub does:
   //   TMatrixT tmp;
   // Doesn't compile here, because we are outside the class?
   // So we'll be explicit:
   TMatrixT<Element1> target;

   return OuterProduct(target,v1,v2);
}

//______________________________________________________________________________
template <class Element1,class Element2,class Element3>
TMatrixT<Element1>
&OuterProduct(TMatrixT<Element1> &target,const TVectorT<Element2> &v1,const TVectorT<Element3> &v2)
{
   // Return the matrix M = v1 * v2'

   target.ResizeTo(v1.GetLwb(), v1.GetUpb(), v2.GetLwb(), v2.GetUpb());

         Element1 *       mp      = target.GetMatrixArray();
   const Element1 * const m_last  = mp + target.GetNoElements();

   const Element2 *       v1p     = v1.GetMatrixArray();
   const Element2 * const v1_last = v1p + v1.GetNrows();

   const Element3 * const v20     = v2.GetMatrixArray();
   const Element3 *       v2p     = v20;
   const Element3 * const v2_last = v2p + v2.GetNrows();

   while (v1p < v1_last) {
      v2p = v20;
      while (v2p < v2_last) {
         *mp++ = *v1p * *v2p++ ;
      }
      v1p++;
  }

  R__ASSERT(v1p == v1_last && mp == m_last && v2p == v2_last);

  return target;
}

//______________________________________________________________________________
template <class Element1, class Element2, class Element3>
Element1 Mult(const TVectorT<Element1> &v1,const TMatrixT<Element2> &m,
              const TVectorT<Element3> &v2)
{
  // Perform v1 * M * v2, a scalar result

   if (gMatrixCheck) {
      if (!AreCompatible(v1, m)) {
         ::Error("Mult", "Vector v1 and matrix m incompatible");
         return 0;
      }
      if (!AreCompatible(m, v2)) {
         ::Error("Mult", "Matrix m and vector v2 incompatible");
         return 0;
      }
   }

   const Element1 *       v1p     = v1.GetMatrixArray();    // first of v1
   const Element1 * const v1_last = v1p + v1.GetNrows();    // last of  v1

   const Element2 *       mp      = m.GetMatrixArray();     // first of m
   const Element2 * const m_last  = mp + m.GetNoElements(); // last of  m

   const Element3 * const v20     = v2.GetMatrixArray();    // first of v2
   const Element3 *       v2p     = v20;                    // running  v2
   const Element3 * const v2_last = v2p + v2.GetNrows();    // last of  v2

   Element1 sum     = 0;      // scalar result accumulator
   Element3 dot     = 0;      // M_row * v2 dot product accumulator

   while (v1p < v1_last) {
      v2p  = v20;               // at beginning of v2
      while (v2p < v2_last) {   // compute (M[i] * v2) dot product
         dot += *mp++ * *v2p++;
      }
      sum += *v1p++ * dot;      // v1[i] * (M[i] * v2)
      dot  = 0;                 // start next dot product
   }

   R__ASSERT(v1p == v1_last && mp == m_last && v2p == v2_last);

   return sum;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &Add(TVectorT<Element> &target,Element scalar,const TVectorT<Element> &source)
{
// Modify addition: target += scalar * source.

   if (gMatrixCheck && !AreCompatible(target,source)) {
      Error("Add(TVectorT<Element> &,Element,const TVectorT<Element> &)","vector's are incompatible");
      return target;
   }

   const Element *       sp  = source.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();
   if (scalar == 1.0 ) {
      while ( tp < ftp )
         *tp++ += *sp++;
   } else if (scalar == -1.0) {
      while ( tp < ftp )
         *tp++ -= *sp++;
   } else {
      while ( tp < ftp )
         *tp++ += scalar * *sp++;
   }

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &Add(TVectorT<Element> &target,Element scalar,
                       const TMatrixT<Element> &a,const TVectorT<Element> &source)
{
// Modify addition: target += scalar * A * source.
// NOTE: in case scalar=0, do  target = A * source.

   if (gMatrixCheck) {
      R__ASSERT(target.IsValid());
      R__ASSERT(a.IsValid());
      if (a.GetNrows() != target.GetNrows() || a.GetRowLwb() != target.GetLwb()) {
         Error("Add(TVectorT &,const TMatrixT &,const TVectorT &)","target vector and matrix are incompatible");
         return target;
      }

      R__ASSERT(source.IsValid());
      if (a.GetNcols() != source.GetNrows() || a.GetColLwb() != source.GetLwb()) {
         Error("Add(TVectorT &,const TMatrixT &,const TVectorT &)","source vector and matrix are incompatible");
         return target;
      }
   }

   const Element * const sp = source.GetMatrixArray();  // sources vector ptr
   const Element *       mp = a.GetMatrixArray();       // Matrix row ptr
         Element *       tp = target.GetMatrixArray();  // Target vector ptr
#ifdef CBLAS
   if (typeid(Element) == typeid(Double_t))
      cblas_dsymv(CblasRowMajor,CblasUpper,fNrows,scalar,mp,
                  fNrows,sp,1,0.0,tp,1);
   else if (typeid(Element) != typeid(Float_t))
      cblas_ssymv(CblasRowMajor,CblasUpper,fNrows,scalar,mp,
                  fNrows,sp,1,0.0,tp,1);
   else
      Error("operator*=","type %s not implemented in BLAS library",typeid(Element));
#else
   const Element * const sp_last = sp+source.GetNrows();
   const Element * const tp_last = tp+target.GetNrows();
   if (scalar == 1.0) {
      while (tp < tp_last) {
         const Element *sp1 = sp;
         Element sum = 0;
         while (sp1 < sp_last)
            sum += *sp1++ * *mp++;
         *tp++ += sum;
      }
   } else if (scalar == 0.0) {
      while (tp < tp_last) {
         const Element *sp1 = sp;
         Element sum = 0;
         while (sp1 < sp_last)
            sum += *sp1++ * *mp++;
         *tp++  = sum;
      }
   } else if (scalar == -1.0) {
      while (tp < tp_last) {
         const Element *sp1 = sp;
         Element sum = 0;
         while (sp1 < sp_last)
            sum += *sp1++ * *mp++;
         *tp++ -= sum;
      }
   } else {
      while (tp < tp_last) {
        const Element *sp1 = sp;
        Element sum = 0;
        while (sp1 < sp_last)
           sum += *sp1++ * *mp++;
        *tp++ += scalar * sum;
      }
   }

   if (gMatrixCheck) R__ASSERT(mp == a.GetMatrixArray()+a.GetNoElements());
#endif

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &Add(TVectorT<Element> &target,Element scalar,
                       const TMatrixTSym<Element> &a,const TVectorT<Element> &source)
{
// Modify addition: target += A * source.
// NOTE: in case scalar=0, do  target = A * source.

   if (gMatrixCheck) {
      R__ASSERT(target.IsValid());
      R__ASSERT(source.IsValid());
      R__ASSERT(a.IsValid());
      if (a.GetNrows() != target.GetNrows() || a.GetRowLwb() != target.GetLwb()) {
         Error("Add(TVectorT &,const TMatrixT &,const TVectorT &)","target vector and matrix are incompatible");
         return target;
      }
   }

   const Element * const sp = source.GetMatrixArray();  // sources vector ptr
   const Element *       mp = a.GetMatrixArray();       // Matrix row ptr
         Element *       tp = target.GetMatrixArray();  // Target vector ptr
#ifdef CBLAS
   if (typeid(Element) == typeid(Double_t))
      cblas_dsymv(CblasRowMajor,CblasUpper,fNrows,1.0,mp,
                  fNrows,sp,1,0.0,tp,1);
   else if (typeid(Element) != typeid(Float_t))
      cblas_ssymv(CblasRowMajor,CblasUpper,fNrows,1.0,mp,
                  fNrows,sp,1,0.0,tp,1);
   else
      Error("operator*=","type %s not implemented in BLAS library",typeid(Element));
#else
   const Element * const sp_last = sp+source.GetNrows();
   const Element * const tp_last = tp+target.GetNrows();
   if (scalar == 1.0) {
      while (tp < tp_last) {
         const Element *sp1 = sp;
         Element sum = 0;
         while (sp1 < sp_last)
            sum += *sp1++ * *mp++;
         *tp++ += sum;
      }
   } else if (scalar == 0.0) {
      while (tp < tp_last) {
         const Element *sp1 = sp;
         Element sum = 0;
         while (sp1 < sp_last)
            sum += *sp1++ * *mp++;
         *tp++  = sum;
      }
   } else if (scalar == -1.0) {
      while (tp < tp_last) {
         const Element *sp1 = sp;
         Element sum = 0;
         while (sp1 < sp_last)
            sum += *sp1++ * *mp++;
         *tp++ -= sum;
      }
   } else {
      while (tp < tp_last) {
         const Element *sp1 = sp;
         Element sum = 0;
         while (sp1 < sp_last)
            sum += *sp1++ * *mp++;
         *tp++ += scalar * sum;
      }
   }
   R__ASSERT(mp == a.GetMatrixArray()+a.GetNoElements());
#endif

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &Add(TVectorT<Element> &target,Element scalar,
                       const TMatrixTSparse<Element> &a,const TVectorT<Element> &source)
{
// Modify addition: target += A * source.
// NOTE: in case scalar=0, do  target = A * source.

   if (gMatrixCheck) {
      R__ASSERT(target.IsValid());
      R__ASSERT(a.IsValid());
      if (a.GetNrows() != target.GetNrows() || a.GetRowLwb() != target.GetLwb()) {
         Error("Add(TVectorT &,const TMatrixT &,const TVectorT &)","target vector and matrix are incompatible");
         return target;
      }

      R__ASSERT(source.IsValid());
      if (a.GetNcols() != source.GetNrows() || a.GetColLwb() != source.GetLwb()) {
         Error("Add(TVectorT &,const TMatrixT &,const TVectorT &)","source vector and matrix are incompatible");
         return target;
      }
   }

   const Int_t   * const pRowIndex = a.GetRowIndexArray();
   const Int_t   * const pColIndex = a.GetColIndexArray();
   const Element * const mp        = a.GetMatrixArray();     // Matrix row ptr

   const Element * const sp = source.GetMatrixArray(); // Source vector ptr
         Element *       tp = target.GetMatrixArray(); // Target vector ptr

   if (scalar == 1.0) {
      for (Int_t irow = 0; irow < a.GetNrows(); irow++) {
         const Int_t sIndex = pRowIndex[irow];
         const Int_t eIndex = pRowIndex[irow+1];
         Element sum = 0.0;
         for (Int_t index = sIndex; index < eIndex; index++) {
            const Int_t icol = pColIndex[index];
            sum += mp[index]*sp[icol];
         }
         tp[irow] += sum;
      }
   } else if (scalar == 0.0) {
      for (Int_t irow = 0; irow < a.GetNrows(); irow++) {
         const Int_t sIndex = pRowIndex[irow];
         const Int_t eIndex = pRowIndex[irow+1];
         Element sum = 0.0;
         for (Int_t index = sIndex; index < eIndex; index++) {
            const Int_t icol = pColIndex[index];
            sum += mp[index]*sp[icol];
         }
         tp[irow]  = sum;
      }
   } else if (scalar == -1.0) {
     for (Int_t irow = 0; irow < a.GetNrows(); irow++) {
        const Int_t sIndex = pRowIndex[irow];
        const Int_t eIndex = pRowIndex[irow+1];
        Element sum = 0.0;
        for (Int_t index = sIndex; index < eIndex; index++) {
           const Int_t icol = pColIndex[index];
           sum += mp[index]*sp[icol];
        }
        tp[irow] -= sum;
      }
   } else {
      for (Int_t irow = 0; irow < a.GetNrows(); irow++) {
        const Int_t sIndex = pRowIndex[irow];
        const Int_t eIndex = pRowIndex[irow+1];
        Element sum = 0.0;
        for (Int_t index = sIndex; index < eIndex; index++) {
           const Int_t icol = pColIndex[index];
           sum += mp[index]*sp[icol];
        }
        tp[irow] += scalar * sum;
      }
   }

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &AddElemMult(TVectorT<Element> &target,Element scalar,
                      const TVectorT<Element> &source1,const TVectorT<Element> &source2)
{
// Modify addition: target += scalar * ElementMult(source1,source2) .

   if (gMatrixCheck && !(AreCompatible(target,source1) && AreCompatible(target,source1))) {
      Error("AddElemMult(TVectorT<Element> &,Element,const TVectorT<Element> &,const TVectorT<Element> &)",
             "vector's are incompatible");
      return target;
   }

   const Element *       sp1 = source1.GetMatrixArray();
   const Element *       sp2 = source2.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();

   if (scalar == 1.0 ) {
      while ( tp < ftp )
         *tp++ += *sp1++ * *sp2++;
   } else if (scalar == -1.0) {
      while ( tp < ftp )
         *tp++ -= *sp1++ * *sp2++;
   } else {
      while ( tp < ftp )
         *tp++ += scalar * *sp1++ * *sp2++;
   }

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &AddElemMult(TVectorT<Element> &target,Element scalar,
                      const TVectorT<Element> &source1,const TVectorT<Element> &source2,const TVectorT<Element> &select)
{
// Modify addition: target += scalar * ElementMult(source1,source2) only for those elements
// where select[i] != 0.0

   if (gMatrixCheck && !( AreCompatible(target,source1) && AreCompatible(target,source1) &&
          AreCompatible(target,select) )) {
      Error("AddElemMult(TVectorT<Element> &,Element,const TVectorT<Element> &,const TVectorT<Element> &,onst TVectorT<Element> &)",
             "vector's are incompatible");
      return target;
   }

   const Element *       sp1 = source1.GetMatrixArray();
   const Element *       sp2 = source2.GetMatrixArray();
   const Element *       mp  = select.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();

   if (scalar == 1.0 ) {
      while ( tp < ftp ) {
         if (*mp) *tp += *sp1 * *sp2;
         mp++; tp++; sp1++; sp2++;
      }
   } else if (scalar == -1.0) {
      while ( tp < ftp ) {
         if (*mp) *tp -= *sp1 * *sp2;
         mp++; tp++; sp1++; sp2++;
      }
   } else {
      while ( tp < ftp ) {
         if (*mp) *tp += scalar * *sp1 * *sp2;
         mp++; tp++; sp1++; sp2++;
      }
   }

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &AddElemDiv(TVectorT<Element> &target,Element scalar,
                     const TVectorT<Element> &source1,const TVectorT<Element> &source2)
{
// Modify addition: target += scalar * ElementDiv(source1,source2) .

   if (gMatrixCheck && !(AreCompatible(target,source1) && AreCompatible(target,source1))) {
      Error("AddElemDiv(TVectorT<Element> &,Element,const TVectorT<Element> &,const TVectorT<Element> &)",
             "vector's are incompatible");
      return target;
   }

   const Element *       sp1 = source1.GetMatrixArray();
   const Element *       sp2 = source2.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();

   if (scalar == 1.0 ) {
      while ( tp < ftp ) {
         if (*sp2 != 0.0)
            *tp += *sp1 / *sp2;
         else {
            const Int_t irow = (sp2-source2.GetMatrixArray())/source2.GetNrows();
            Error("AddElemDiv","source2 (%d) is zero",irow);
         }
         tp++; sp1++; sp2++;
      }
   } else if (scalar == -1.0) {
      while ( tp < ftp ) {
         if (*sp2 != 0.0)
            *tp -= *sp1 / *sp2;
         else {
            const Int_t irow = (sp2-source2.GetMatrixArray())/source2.GetNrows();
            Error("AddElemDiv","source2 (%d) is zero",irow);
         }
         tp++; sp1++; sp2++;
      }
   } else {
      while ( tp < ftp ) {
         if (*sp2 != 0.0)
            *tp += scalar * *sp1 / *sp2;
         else {
            const Int_t irow = (sp2-source2.GetMatrixArray())/source2.GetNrows();
            Error("AddElemDiv","source2 (%d) is zero",irow);
         }
         tp++; sp1++; sp2++;
      }
   }

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &AddElemDiv(TVectorT<Element> &target,Element scalar,
                     const TVectorT<Element> &source1,const TVectorT<Element> &source2,const TVectorT<Element> &select)
{
// Modify addition: target += scalar * ElementDiv(source1,source2) only for those elements
// where select[i] != 0.0

   if (gMatrixCheck && !( AreCompatible(target,source1) && AreCompatible(target,source1) &&
          AreCompatible(target,select) )) {
      Error("AddElemDiv(TVectorT<Element> &,Element,const TVectorT<Element> &,const TVectorT<Element> &,onst TVectorT<Element> &)",
             "vector's are incompatible");
      return target;
   }

   const Element *       sp1 = source1.GetMatrixArray();
   const Element *       sp2 = source2.GetMatrixArray();
   const Element *       mp  = select.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();

   if (scalar == 1.0 ) {
      while ( tp < ftp ) {
         if (*mp) {
            if (*sp2 != 0.0)
               *tp += *sp1 / *sp2;
            else {
               const Int_t irow = (sp2-source2.GetMatrixArray())/source2.GetNrows();
               Error("AddElemDiv","source2 (%d) is zero",irow);
            }
         }
         mp++; tp++; sp1++; sp2++;
      }
   } else if (scalar == -1.0) {
      while ( tp < ftp ) {
         if (*mp) {
            if (*sp2 != 0.0)
               *tp -= *sp1 / *sp2;
            else {
               const Int_t irow = (sp2-source2.GetMatrixArray())/source2.GetNrows();
               Error("AddElemDiv","source2 (%d) is zero",irow);
            }
         }
         mp++; tp++; sp1++; sp2++;
      }
   } else {
      while ( tp < ftp ) {
         if (*mp) {
            if (*sp2 != 0.0)
               *tp += scalar * *sp1 / *sp2;
            else {
               const Int_t irow = (sp2-source2.GetMatrixArray())/source2.GetNrows();
               Error("AddElemDiv","source2 (%d) is zero",irow);
            }
         }
         mp++; tp++; sp1++; sp2++;
      }
   }

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &ElementMult(TVectorT<Element> &target,const TVectorT<Element> &source)
{
// Multiply target by the source, element-by-element.

   if (gMatrixCheck && !AreCompatible(target,source)) {
      Error("ElementMult(TVectorT<Element> &,const TVectorT<Element> &)","vector's are incompatible");
      return target;
   }

   const Element *       sp  = source.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();
   while ( tp < ftp )
      *tp++ *= *sp++;

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &ElementMult(TVectorT<Element> &target,const TVectorT<Element> &source,const TVectorT<Element> &select)
{
// Multiply target by the source, element-by-element only where select[i] != 0.0

   if (gMatrixCheck && !(AreCompatible(target,source) && AreCompatible(target,select))) {
      Error("ElementMult(TVectorT<Element> &,const TVectorT<Element> &,const TVectorT<Element> &)","vector's are incompatible");
      return target;
   }

   const Element *       sp  = source.GetMatrixArray();
   const Element *       mp  = select.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();
   while ( tp < ftp ) {
      if (*mp) *tp *= *sp;
      mp++; tp++; sp++;
   }

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &ElementDiv(TVectorT<Element> &target,const TVectorT<Element> &source)
{
// Divide target by the source, element-by-element.

   if (gMatrixCheck && !AreCompatible(target,source)) {
      Error("ElementDiv(TVectorT<Element> &,const TVectorT<Element> &)","vector's are incompatible");
      return target;
   }

   const Element *       sp  = source.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();
   while ( tp < ftp ) {
      if (*sp  != 0.0)
         *tp++ /= *sp++;
      else {
         const Int_t irow = (sp-source.GetMatrixArray())/source.GetNrows();
         Error("ElementDiv","source (%d) is zero",irow);
      }
   }

   return target;
}

//______________________________________________________________________________
template<class Element>
TVectorT<Element> &ElementDiv(TVectorT<Element> &target,const TVectorT<Element> &source,const TVectorT<Element> &select)
{
// Divide target by the source, element-by-element only where select[i] != 0.0

   if (gMatrixCheck && !AreCompatible(target,source)) {
      Error("ElementDiv(TVectorT<Element> &,const TVectorT<Element> &,const TVectorT<Element> &)","vector's are incompatible");
      return target;
   }

   const Element *       sp  = source.GetMatrixArray();
   const Element *       mp  = select.GetMatrixArray();
         Element *       tp  = target.GetMatrixArray();
   const Element * const ftp = tp+target.GetNrows();
   while ( tp < ftp ) {
      if (*mp) {
         if (*sp != 0.0)
            *tp /= *sp;
         else {
            const Int_t irow = (sp-source.GetMatrixArray())/source.GetNrows();
            Error("ElementDiv","source (%d) is zero",irow);
         }
      }
      mp++; tp++; sp++;
   }

   return target;
}

//______________________________________________________________________________
template<class Element1,class Element2>
Bool_t AreCompatible(const TVectorT<Element1> &v1,const TVectorT<Element2> &v2,Int_t verbose)
{
// Check if v1 and v2 are both valid and have the same shape

   if (!v1.IsValid()) {
      if (verbose)
         ::Error("AreCompatible", "vector 1 not valid");
      return kFALSE;
   }
   if (!v2.IsValid()) {
      if (verbose)
         ::Error("AreCompatible", "vector 2 not valid");
      return kFALSE;
   }

   if (v1.GetNrows() != v2.GetNrows() || v1.GetLwb() != v2.GetLwb()) {
      if (verbose)
         ::Error("AreCompatible", "matrices 1 and 2 not compatible");
      return kFALSE;
   }

   return kTRUE;
}

//______________________________________________________________________________
template<class Element1, class Element2>
Bool_t AreCompatible(const TMatrixT<Element1> &m,const TVectorT<Element2> &v,Int_t verbose)
{
  // Check if m and v are both valid and have compatible shapes for M * v

   if (!m.IsValid()) {
      if (verbose)
         ::Error("AreCompatible", "Matrix not valid");
      return kFALSE;
   }
   if (!v.IsValid()) {
      if (verbose)
         ::Error("AreCompatible", "vector not valid");
      return kFALSE;
   }

   if (m.GetNcols() != v.GetNrows() ) {
      if (verbose)
         ::Error("AreCompatible", "matrix and vector not compatible");
      return kFALSE;
   }

   return kTRUE;
}

//______________________________________________________________________________
template<class Element1, class Element2>
Bool_t AreCompatible(const TVectorT<Element1> &v,const TMatrixT<Element2> &m,Int_t verbose)
{
  // Check if m and v are both valid and have compatible shapes for v * M

   if (!m.IsValid()) {
      if (verbose)
         ::Error("AreCompatible", "Matrix not valid");
      return kFALSE;
   }
   if (!v.IsValid()) {
      if (verbose)
         ::Error("AreCompatible", "vector not valid");
      return kFALSE;
   }

   if (v.GetNrows() != m.GetNrows() ) {
      if (verbose)
         ::Error("AreCompatible", "vector and matrix not compatible");
      return kFALSE;
   }

   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
void Compare(const TVectorT<Element> &v1,const TVectorT<Element> &v2)
{
// Compare two vectors and print out the result of the comparison.

   if (!AreCompatible(v1,v2)) {
      Error("Compare(const TVectorT<Element> &,const TVectorT<Element> &)","vectors are incompatible");
      return;
   }

   printf("\n\nComparison of two TVectorTs:\n");

   Element norm1  = 0;       // Norm of the Matrices
   Element norm2  = 0;       // Norm of the Matrices
   Element ndiff  = 0;       // Norm of the difference
   Int_t    imax   = 0;       // For the elements that differ most
   Element difmax = -1;
   const Element *mp1 = v1.GetMatrixArray();    // Vector element pointers
   const Element *mp2 = v2.GetMatrixArray();

   for (Int_t i = 0; i < v1.GetNrows(); i++) {
      const Element mv1  = *mp1++;
      const Element mv2  = *mp2++;
      const Element diff = TMath::Abs(mv1-mv2);

      if (diff > difmax) {
         difmax = diff;
         imax = i;
      }
      norm1 += TMath::Abs(mv1);
      norm2 += TMath::Abs(mv2);
      ndiff += TMath::Abs(diff);
   }

   imax += v1.GetLwb();
   printf("\nMaximal discrepancy    \t\t%g",difmax);
   printf("\n   occured at the point\t\t(%d)",imax);
   const Element mv1 = v1(imax);
   const Element mv2 = v2(imax);
   printf("\n Vector 1 element is    \t\t%g",mv1);
   printf("\n Vector 2 element is    \t\t%g",mv2);
   printf("\n Absolute error v2[i]-v1[i]\t\t%g",mv2-mv1);
   printf("\n Relative error\t\t\t\t%g\n",
          (mv2-mv1)/TMath::Max(TMath::Abs(mv2+mv1)/2,(Element)1e-7));

   printf("\n||Vector 1||   \t\t\t%g",norm1);
   printf("\n||Vector 2||   \t\t\t%g",norm2);
   printf("\n||Vector1-Vector2||\t\t\t\t%g",ndiff);
   printf("\n||Vector1-Vector2||/sqrt(||Vector1|| ||Vector2||)\t%g\n\n",
          ndiff/TMath::Max(TMath::Sqrt(norm1*norm2),1e-7));
}

//______________________________________________________________________________
template<class Element>
Bool_t VerifyVectorValue(const TVectorT<Element> &v,Element val,
                         Int_t verbose,Element maxDevAllow)
{
// Validate that all elements of vector have value val within maxDevAllow .

   Int_t   imax      = 0;
   Element maxDevObs = 0;

   if (TMath::Abs(maxDevAllow) <= 0.0)
      maxDevAllow = std::numeric_limits<Element>::epsilon();

   for (Int_t i = v.GetLwb(); i <= v.GetUpb(); i++) {
      const Element dev = TMath::Abs(v(i)-val);
      if (dev > maxDevObs) {
         imax      = i;
         maxDevObs = dev;
      }
   }

   if (maxDevObs == 0)
      return kTRUE;

   if (verbose) {
      printf("Largest dev for (%d); dev = |%g - %g| = %g\n",imax,v(imax),val,maxDevObs);
      if (maxDevObs > maxDevAllow)
         Error("VerifyVectorValue","Deviation > %g\n",maxDevAllow);
   }

   if (maxDevObs > maxDevAllow)
      return kFALSE;
   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
Bool_t VerifyVectorIdentity(const TVectorT<Element> &v1,const TVectorT<Element> &v2,
                            Int_t verbose, Element maxDevAllow)
{
// Verify that elements of the two vectors are equal within maxDevAllow .

   Int_t   imax      = 0;
   Element maxDevObs = 0;

   if (!AreCompatible(v1,v2))
      return kFALSE;

   if (TMath::Abs(maxDevAllow) <= 0.0)
      maxDevAllow = std::numeric_limits<Element>::epsilon();

   for (Int_t i = v1.GetLwb(); i <= v1.GetUpb(); i++) {
      const Element dev = TMath::Abs(v1(i)-v2(i));
      if (dev > maxDevObs) {
         imax      = i;
         maxDevObs = dev;
      }
   }

   if (maxDevObs == 0)
      return kTRUE;

   if (verbose) {
      printf("Largest dev for (%d); dev = |%g - %g| = %g\n",imax,v1(imax),v2(imax),maxDevObs);
      if (maxDevObs > maxDevAllow)
         Error("VerifyVectorIdentity","Deviation > %g\n",maxDevAllow);
   }

   if (maxDevObs > maxDevAllow) {
      return kFALSE;
   }
   return kTRUE;
}

//______________________________________________________________________________
template<class Element>
void TVectorT<Element>::Streamer(TBuffer &R__b)
{
// Stream an object of class TVectorT.

   if (R__b.IsReading()) {
      UInt_t R__s, R__c;
      Version_t R__v = R__b.ReadVersion(&R__s,&R__c);
      if (R__v > 1) {
         Clear();
         R__b.ReadClassBuffer(TVectorT<Element>::Class(),this,R__v,R__s,R__c);
      } else { //====process old versions before automatic schema evolution
         TObject::Streamer(R__b);
         R__b >> fRowLwb;
         fNrows = R__b.ReadArray(fElements);
         R__b.CheckByteCount(R__s, R__c, TVectorT<Element>::IsA());
      }
      if (fNrows > 0 && fNrows <= kSizeMax) {
         memcpy(fDataStack,fElements,fNrows*sizeof(Element));
         delete [] fElements;
         fElements = fDataStack;
      } else if (fNrows < 0)
         Invalidate();

      if (R__v < 3)
         MakeValid();
    } else {
       R__b.WriteClassBuffer(TVectorT<Element>::Class(),this);
   }
}

#ifndef ROOT_TMatrixFfwd
#include "TMatrixFfwd.h"
#endif
#ifndef ROOT_TMatrixFSymfwd
#include "TMatrixFSymfwd.h"
#endif
#ifndef ROOT_TMatrixFSparsefwd
#include "TMatrixFSparsefwd.h"
#endif

template class TVectorT<Float_t>;

template Bool_t    operator==          <Float_t>(const TVectorF       &source1,const TVectorF &source2);
template TVectorF  operator+           <Float_t>(const TVectorF       &source1,const TVectorF &source2);
template TVectorF  operator-           <Float_t>(const TVectorF       &source1,const TVectorF &source2);
template Float_t   operator*           <Float_t>(const TVectorF       &source1,const TVectorF &source2);
template TVectorF  operator*           <Float_t>(const TMatrixF       &a,      const TVectorF &source);
template TVectorF  operator*           <Float_t>(const TMatrixFSym    &a,      const TVectorF &source);
template TVectorF  operator*           <Float_t>(const TMatrixFSparse &a,      const TVectorF &source);
template TVectorF  operator*           <Float_t>(      Float_t         val,    const TVectorF &source);

template Float_t   Dot                 <Float_t>(const TVectorF       &v1,     const TVectorF &v2);
template TMatrixF  OuterProduct        <Float_t,Float_t>
                                                (const TVectorF       &v1,     const TVectorF &v2);
template TMatrixF &OuterProduct        <Float_t,Float_t,Float_t>
                                                (      TMatrixF       &target, const TVectorF &v1,     const TVectorF       &v2);
template Float_t   Mult                <Float_t,Float_t,Float_t>
                                                (const TVectorF       &v1,     const TMatrixF &m,      const TVectorF       &v2);

template TVectorF &Add                 <Float_t>(      TVectorF       &target,       Float_t   scalar, const TVectorF       &source);
template TVectorF &Add                 <Float_t>(      TVectorF       &target,       Float_t   scalar, const TMatrixF       &a,
                                                                               const TVectorF &source);
template TVectorF &Add                 <Float_t>(      TVectorF       &target,       Float_t   scalar, const TMatrixFSym    &a,
                                                                               const TVectorF &source);
template TVectorF &Add                 <Float_t>(      TVectorF       &target,       Float_t   scalar, const TMatrixFSparse &a,
                                                                               const TVectorF &source);
template TVectorF &AddElemMult         <Float_t>(      TVectorF       &target,       Float_t   scalar, const TVectorF       &source1,
                                                                               const TVectorF &source2);
template TVectorF &AddElemMult         <Float_t>(      TVectorF       &target,       Float_t   scalar, const TVectorF       &source1,
                                                                               const TVectorF &source2,const TVectorF       &select);
template TVectorF &AddElemDiv          <Float_t>(      TVectorF       &target,       Float_t   scalar, const TVectorF       &source1,
                                                                               const TVectorF &source2);
template TVectorF &AddElemDiv          <Float_t>(      TVectorF       &target,       Float_t   scalar, const TVectorF       &source1,
                                                                               const TVectorF &source2,const TVectorF       &select);
template TVectorF &ElementMult         <Float_t>(      TVectorF       &target, const TVectorF &source);
template TVectorF &ElementMult         <Float_t>(      TVectorF       &target, const TVectorF &source, const TVectorF       &select);
template TVectorF &ElementDiv          <Float_t>(      TVectorF       &target, const TVectorF &source);
template TVectorF &ElementDiv          <Float_t>(      TVectorF       &target, const TVectorF &source, const TVectorF       &select);

template Bool_t    AreCompatible       <Float_t,Float_t> (const TVectorF &v1,const TVectorF &v2,Int_t verbose);
template Bool_t    AreCompatible       <Float_t,Double_t>(const TVectorF &v1,const TVectorD &v2,Int_t verbose);
template Bool_t    AreCompatible       <Float_t,Float_t> (const TMatrixF &m, const TVectorF &v, Int_t verbose);
template Bool_t    AreCompatible       <Float_t,Float_t> (const TVectorF &v, const TMatrixF &m, Int_t verbose);

template void      Compare             <Float_t>         (const TVectorF &v1,const TVectorF &v2);
template Bool_t    VerifyVectorValue   <Float_t>         (const TVectorF &m,       Float_t   val,Int_t verbose,Float_t maxDevAllow);
template Bool_t    VerifyVectorIdentity<Float_t>         (const TVectorF &m1,const TVectorF &m2, Int_t verbose,Float_t maxDevAllow);

#ifndef ROOT_TMatrixDfwd
#include "TMatrixDfwd.h"
#endif
#ifndef ROOT_TMatrixDSymfwd
#include "TMatrixDSymfwd.h"
#endif
#ifndef ROOT_TMatrixDSparsefwd
#include "TMatrixDSparsefwd.h"
#endif

template class TVectorT<Double_t>;

template Bool_t    operator==          <Double_t>(const TVectorD       &source1,const TVectorD &source2);
template TVectorD  operator+           <Double_t>(const TVectorD       &source1,const TVectorD &source2);
template TVectorD  operator-           <Double_t>(const TVectorD       &source1,const TVectorD &source2);
template Double_t  operator*           <Double_t>(const TVectorD       &source1,const TVectorD &source2);
template TVectorD  operator*           <Double_t>(const TMatrixD       &a,      const TVectorD &source);
template TVectorD  operator*           <Double_t>(const TMatrixDSym    &a,      const TVectorD &source);
template TVectorD  operator*           <Double_t>(const TMatrixDSparse &a,      const TVectorD &source);
template TVectorD  operator*           <Double_t>(      Double_t        val,    const TVectorD &source);

template Double_t  Dot                 <Double_t>(const TVectorD       &v1,     const TVectorD &v2);
template TMatrixD  OuterProduct        <Double_t,Double_t>
                                                 (const TVectorD       &v1,     const TVectorD &v2);
template TMatrixD &OuterProduct        <Double_t,Double_t,Double_t>
                                                 (      TMatrixD       &target, const TVectorD &v1,     const TVectorD       &v2);
template Double_t  Mult                <Double_t,Double_t,Double_t>
                                                 (const TVectorD       &v1,     const TMatrixD &m,      const TVectorD       &v2);

template TVectorD &Add                 <Double_t>(      TVectorD       &target,       Double_t  scalar, const TVectorD       &source);
template TVectorD &Add                 <Double_t>(      TVectorD       &target,       Double_t  scalar, const TMatrixD       &a,
                                                                                const TVectorD &source);
template TVectorD &Add                 <Double_t>(      TVectorD       &target,       Double_t  scalar, const TMatrixDSym    &a
                                                                                ,      const TVectorD &source);
template TVectorD &Add                 <Double_t>(      TVectorD       &target,       Double_t  scalar, const TMatrixDSparse &a
                                                                                ,      const TVectorD &source);
template TVectorD &AddElemMult         <Double_t>(      TVectorD       &target,       Double_t  scalar, const TVectorD       &source1,
                                                                                const TVectorD &source2);
template TVectorD &AddElemMult         <Double_t>(      TVectorD       &target,       Double_t  scalar, const TVectorD       &source1,
                                                                                const TVectorD &source2,const TVectorD       &select);
template TVectorD &AddElemDiv          <Double_t>(      TVectorD       &target,       Double_t  scalar, const TVectorD       &source1,
                                                                                const TVectorD &source2);
template TVectorD &AddElemDiv          <Double_t>(      TVectorD       &target,       Double_t  scalar, const TVectorD       &source1,
                                                                                const TVectorD &source2,const TVectorD       &select);
template TVectorD &ElementMult         <Double_t>(      TVectorD       &target, const TVectorD &source);
template TVectorD &ElementMult         <Double_t>(      TVectorD       &target, const TVectorD &source, const TVectorD       &select);
template TVectorD &ElementDiv          <Double_t>(      TVectorD       &target, const TVectorD &source);
template TVectorD &ElementDiv          <Double_t>(      TVectorD       &target, const TVectorD &source, const TVectorD       &select);

template Bool_t    AreCompatible       <Double_t,Double_t>(const TVectorD &v1,const TVectorD &v2,Int_t verbose);
template Bool_t    AreCompatible       <Double_t,Float_t> (const TVectorD &v1,const TVectorF &v2,Int_t verbose);
template Bool_t    AreCompatible       <Double_t,Double_t>(const TMatrixD &m, const TVectorD &v, Int_t verbose);
template Bool_t    AreCompatible       <Double_t,Double_t>(const TVectorD &v, const TMatrixD &m, Int_t verbose);

template void      Compare             <Double_t>         (const TVectorD &v1,const TVectorD &v2);
template Bool_t    VerifyVectorValue   <Double_t>         (const TVectorD &m,       Double_t  val,Int_t verbose,Double_t maxDevAllow);
template Bool_t    VerifyVectorIdentity<Double_t>         (const TVectorD &m1,const TVectorD &m2, Int_t verbose,Double_t maxDevAllow);
 TVectorT.cxx:1
 TVectorT.cxx:2
 TVectorT.cxx:3
 TVectorT.cxx:4
 TVectorT.cxx:5
 TVectorT.cxx:6
 TVectorT.cxx:7
 TVectorT.cxx:8
 TVectorT.cxx:9
 TVectorT.cxx:10
 TVectorT.cxx:11
 TVectorT.cxx:12
 TVectorT.cxx:13
 TVectorT.cxx:14
 TVectorT.cxx:15
 TVectorT.cxx:16
 TVectorT.cxx:17
 TVectorT.cxx:18
 TVectorT.cxx:19
 TVectorT.cxx:20
 TVectorT.cxx:21
 TVectorT.cxx:22
 TVectorT.cxx:23
 TVectorT.cxx:24
 TVectorT.cxx:25
 TVectorT.cxx:26
 TVectorT.cxx:27
 TVectorT.cxx:28
 TVectorT.cxx:29
 TVectorT.cxx:30
 TVectorT.cxx:31
 TVectorT.cxx:32
 TVectorT.cxx:33
 TVectorT.cxx:34
 TVectorT.cxx:35
 TVectorT.cxx:36
 TVectorT.cxx:37
 TVectorT.cxx:38
 TVectorT.cxx:39
 TVectorT.cxx:40
 TVectorT.cxx:41
 TVectorT.cxx:42
 TVectorT.cxx:43
 TVectorT.cxx:44
 TVectorT.cxx:45
 TVectorT.cxx:46
 TVectorT.cxx:47
 TVectorT.cxx:48
 TVectorT.cxx:49
 TVectorT.cxx:50
 TVectorT.cxx:51
 TVectorT.cxx:52
 TVectorT.cxx:53
 TVectorT.cxx:54
 TVectorT.cxx:55
 TVectorT.cxx:56
 TVectorT.cxx:57
 TVectorT.cxx:58
 TVectorT.cxx:59
 TVectorT.cxx:60
 TVectorT.cxx:61
 TVectorT.cxx:62
 TVectorT.cxx:63
 TVectorT.cxx:64
 TVectorT.cxx:65
 TVectorT.cxx:66
 TVectorT.cxx:67
 TVectorT.cxx:68
 TVectorT.cxx:69
 TVectorT.cxx:70
 TVectorT.cxx:71
 TVectorT.cxx:72
 TVectorT.cxx:73
 TVectorT.cxx:74
 TVectorT.cxx:75
 TVectorT.cxx:76
 TVectorT.cxx:77
 TVectorT.cxx:78
 TVectorT.cxx:79
 TVectorT.cxx:80
 TVectorT.cxx:81
 TVectorT.cxx:82
 TVectorT.cxx:83
 TVectorT.cxx:84
 TVectorT.cxx:85
 TVectorT.cxx:86
 TVectorT.cxx:87
 TVectorT.cxx:88
 TVectorT.cxx:89
 TVectorT.cxx:90
 TVectorT.cxx:91
 TVectorT.cxx:92
 TVectorT.cxx:93
 TVectorT.cxx:94
 TVectorT.cxx:95
 TVectorT.cxx:96
 TVectorT.cxx:97
 TVectorT.cxx:98
 TVectorT.cxx:99
 TVectorT.cxx:100
 TVectorT.cxx:101
 TVectorT.cxx:102
 TVectorT.cxx:103
 TVectorT.cxx:104
 TVectorT.cxx:105
 TVectorT.cxx:106
 TVectorT.cxx:107
 TVectorT.cxx:108
 TVectorT.cxx:109
 TVectorT.cxx:110
 TVectorT.cxx:111
 TVectorT.cxx:112
 TVectorT.cxx:113
 TVectorT.cxx:114
 TVectorT.cxx:115
 TVectorT.cxx:116
 TVectorT.cxx:117
 TVectorT.cxx:118
 TVectorT.cxx:119
 TVectorT.cxx:120
 TVectorT.cxx:121
 TVectorT.cxx:122
 TVectorT.cxx:123
 TVectorT.cxx:124
 TVectorT.cxx:125
 TVectorT.cxx:126
 TVectorT.cxx:127
 TVectorT.cxx:128
 TVectorT.cxx:129
 TVectorT.cxx:130
 TVectorT.cxx:131
 TVectorT.cxx:132
 TVectorT.cxx:133
 TVectorT.cxx:134
 TVectorT.cxx:135
 TVectorT.cxx:136
 TVectorT.cxx:137
 TVectorT.cxx:138
 TVectorT.cxx:139
 TVectorT.cxx:140
 TVectorT.cxx:141
 TVectorT.cxx:142
 TVectorT.cxx:143
 TVectorT.cxx:144
 TVectorT.cxx:145
 TVectorT.cxx:146
 TVectorT.cxx:147
 TVectorT.cxx:148
 TVectorT.cxx:149
 TVectorT.cxx:150
 TVectorT.cxx:151
 TVectorT.cxx:152
 TVectorT.cxx:153
 TVectorT.cxx:154
 TVectorT.cxx:155
 TVectorT.cxx:156
 TVectorT.cxx:157
 TVectorT.cxx:158
 TVectorT.cxx:159
 TVectorT.cxx:160
 TVectorT.cxx:161
 TVectorT.cxx:162
 TVectorT.cxx:163
 TVectorT.cxx:164
 TVectorT.cxx:165
 TVectorT.cxx:166
 TVectorT.cxx:167
 TVectorT.cxx:168
 TVectorT.cxx:169
 TVectorT.cxx:170
 TVectorT.cxx:171
 TVectorT.cxx:172
 TVectorT.cxx:173
 TVectorT.cxx:174
 TVectorT.cxx:175
 TVectorT.cxx:176
 TVectorT.cxx:177
 TVectorT.cxx:178
 TVectorT.cxx:179
 TVectorT.cxx:180
 TVectorT.cxx:181
 TVectorT.cxx:182
 TVectorT.cxx:183
 TVectorT.cxx:184
 TVectorT.cxx:185
 TVectorT.cxx:186
 TVectorT.cxx:187
 TVectorT.cxx:188
 TVectorT.cxx:189
 TVectorT.cxx:190
 TVectorT.cxx:191
 TVectorT.cxx:192
 TVectorT.cxx:193
 TVectorT.cxx:194
 TVectorT.cxx:195
 TVectorT.cxx:196
 TVectorT.cxx:197
 TVectorT.cxx:198
 TVectorT.cxx:199
 TVectorT.cxx:200
 TVectorT.cxx:201
 TVectorT.cxx:202
 TVectorT.cxx:203
 TVectorT.cxx:204
 TVectorT.cxx:205
 TVectorT.cxx:206
 TVectorT.cxx:207
 TVectorT.cxx:208
 TVectorT.cxx:209
 TVectorT.cxx:210
 TVectorT.cxx:211
 TVectorT.cxx:212
 TVectorT.cxx:213
 TVectorT.cxx:214
 TVectorT.cxx:215
 TVectorT.cxx:216
 TVectorT.cxx:217
 TVectorT.cxx:218
 TVectorT.cxx:219
 TVectorT.cxx:220
 TVectorT.cxx:221
 TVectorT.cxx:222
 TVectorT.cxx:223
 TVectorT.cxx:224
 TVectorT.cxx:225
 TVectorT.cxx:226
 TVectorT.cxx:227
 TVectorT.cxx:228
 TVectorT.cxx:229
 TVectorT.cxx:230
 TVectorT.cxx:231
 TVectorT.cxx:232
 TVectorT.cxx:233
 TVectorT.cxx:234
 TVectorT.cxx:235
 TVectorT.cxx:236
 TVectorT.cxx:237
 TVectorT.cxx:238
 TVectorT.cxx:239
 TVectorT.cxx:240
 TVectorT.cxx:241
 TVectorT.cxx:242
 TVectorT.cxx:243
 TVectorT.cxx:244
 TVectorT.cxx:245
 TVectorT.cxx:246
 TVectorT.cxx:247
 TVectorT.cxx:248
 TVectorT.cxx:249
 TVectorT.cxx:250
 TVectorT.cxx:251
 TVectorT.cxx:252
 TVectorT.cxx:253
 TVectorT.cxx:254
 TVectorT.cxx:255
 TVectorT.cxx:256
 TVectorT.cxx:257
 TVectorT.cxx:258
 TVectorT.cxx:259
 TVectorT.cxx:260
 TVectorT.cxx:261
 TVectorT.cxx:262
 TVectorT.cxx:263
 TVectorT.cxx:264
 TVectorT.cxx:265
 TVectorT.cxx:266
 TVectorT.cxx:267
 TVectorT.cxx:268
 TVectorT.cxx:269
 TVectorT.cxx:270
 TVectorT.cxx:271
 TVectorT.cxx:272
 TVectorT.cxx:273
 TVectorT.cxx:274
 TVectorT.cxx:275
 TVectorT.cxx:276
 TVectorT.cxx:277
 TVectorT.cxx:278
 TVectorT.cxx:279
 TVectorT.cxx:280
 TVectorT.cxx:281
 TVectorT.cxx:282
 TVectorT.cxx:283
 TVectorT.cxx:284
 TVectorT.cxx:285
 TVectorT.cxx:286
 TVectorT.cxx:287
 TVectorT.cxx:288
 TVectorT.cxx:289
 TVectorT.cxx:290
 TVectorT.cxx:291
 TVectorT.cxx:292
 TVectorT.cxx:293
 TVectorT.cxx:294
 TVectorT.cxx:295
 TVectorT.cxx:296
 TVectorT.cxx:297
 TVectorT.cxx:298
 TVectorT.cxx:299
 TVectorT.cxx:300
 TVectorT.cxx:301
 TVectorT.cxx:302
 TVectorT.cxx:303
 TVectorT.cxx:304
 TVectorT.cxx:305
 TVectorT.cxx:306
 TVectorT.cxx:307
 TVectorT.cxx:308
 TVectorT.cxx:309
 TVectorT.cxx:310
 TVectorT.cxx:311
 TVectorT.cxx:312
 TVectorT.cxx:313
 TVectorT.cxx:314
 TVectorT.cxx:315
 TVectorT.cxx:316
 TVectorT.cxx:317
 TVectorT.cxx:318
 TVectorT.cxx:319
 TVectorT.cxx:320
 TVectorT.cxx:321
 TVectorT.cxx:322
 TVectorT.cxx:323
 TVectorT.cxx:324
 TVectorT.cxx:325
 TVectorT.cxx:326
 TVectorT.cxx:327
 TVectorT.cxx:328
 TVectorT.cxx:329
 TVectorT.cxx:330
 TVectorT.cxx:331
 TVectorT.cxx:332
 TVectorT.cxx:333
 TVectorT.cxx:334
 TVectorT.cxx:335
 TVectorT.cxx:336
 TVectorT.cxx:337
 TVectorT.cxx:338
 TVectorT.cxx:339
 TVectorT.cxx:340
 TVectorT.cxx:341
 TVectorT.cxx:342
 TVectorT.cxx:343
 TVectorT.cxx:344
 TVectorT.cxx:345
 TVectorT.cxx:346
 TVectorT.cxx:347
 TVectorT.cxx:348
 TVectorT.cxx:349
 TVectorT.cxx:350
 TVectorT.cxx:351
 TVectorT.cxx:352
 TVectorT.cxx:353
 TVectorT.cxx:354
 TVectorT.cxx:355
 TVectorT.cxx:356
 TVectorT.cxx:357
 TVectorT.cxx:358
 TVectorT.cxx:359
 TVectorT.cxx:360
 TVectorT.cxx:361
 TVectorT.cxx:362
 TVectorT.cxx:363
 TVectorT.cxx:364
 TVectorT.cxx:365
 TVectorT.cxx:366
 TVectorT.cxx:367
 TVectorT.cxx:368
 TVectorT.cxx:369
 TVectorT.cxx:370
 TVectorT.cxx:371
 TVectorT.cxx:372
 TVectorT.cxx:373
 TVectorT.cxx:374
 TVectorT.cxx:375
 TVectorT.cxx:376
 TVectorT.cxx:377
 TVectorT.cxx:378
 TVectorT.cxx:379
 TVectorT.cxx:380
 TVectorT.cxx:381
 TVectorT.cxx:382
 TVectorT.cxx:383
 TVectorT.cxx:384
 TVectorT.cxx:385
 TVectorT.cxx:386
 TVectorT.cxx:387
 TVectorT.cxx:388
 TVectorT.cxx:389
 TVectorT.cxx:390
 TVectorT.cxx:391
 TVectorT.cxx:392
 TVectorT.cxx:393
 TVectorT.cxx:394
 TVectorT.cxx:395
 TVectorT.cxx:396
 TVectorT.cxx:397
 TVectorT.cxx:398
 TVectorT.cxx:399
 TVectorT.cxx:400
 TVectorT.cxx:401
 TVectorT.cxx:402
 TVectorT.cxx:403
 TVectorT.cxx:404
 TVectorT.cxx:405
 TVectorT.cxx:406
 TVectorT.cxx:407
 TVectorT.cxx:408
 TVectorT.cxx:409
 TVectorT.cxx:410
 TVectorT.cxx:411
 TVectorT.cxx:412
 TVectorT.cxx:413
 TVectorT.cxx:414
 TVectorT.cxx:415
 TVectorT.cxx:416
 TVectorT.cxx:417
 TVectorT.cxx:418
 TVectorT.cxx:419
 TVectorT.cxx:420
 TVectorT.cxx:421
 TVectorT.cxx:422
 TVectorT.cxx:423
 TVectorT.cxx:424
 TVectorT.cxx:425
 TVectorT.cxx:426
 TVectorT.cxx:427
 TVectorT.cxx:428
 TVectorT.cxx:429
 TVectorT.cxx:430
 TVectorT.cxx:431
 TVectorT.cxx:432
 TVectorT.cxx:433
 TVectorT.cxx:434
 TVectorT.cxx:435
 TVectorT.cxx:436
 TVectorT.cxx:437
 TVectorT.cxx:438
 TVectorT.cxx:439
 TVectorT.cxx:440
 TVectorT.cxx:441
 TVectorT.cxx:442
 TVectorT.cxx:443
 TVectorT.cxx:444
 TVectorT.cxx:445
 TVectorT.cxx:446
 TVectorT.cxx:447
 TVectorT.cxx:448
 TVectorT.cxx:449
 TVectorT.cxx:450
 TVectorT.cxx:451
 TVectorT.cxx:452
 TVectorT.cxx:453
 TVectorT.cxx:454
 TVectorT.cxx:455
 TVectorT.cxx:456
 TVectorT.cxx:457
 TVectorT.cxx:458
 TVectorT.cxx:459
 TVectorT.cxx:460
 TVectorT.cxx:461
 TVectorT.cxx:462
 TVectorT.cxx:463
 TVectorT.cxx:464
 TVectorT.cxx:465
 TVectorT.cxx:466
 TVectorT.cxx:467
 TVectorT.cxx:468
 TVectorT.cxx:469
 TVectorT.cxx:470
 TVectorT.cxx:471
 TVectorT.cxx:472
 TVectorT.cxx:473
 TVectorT.cxx:474
 TVectorT.cxx:475
 TVectorT.cxx:476
 TVectorT.cxx:477
 TVectorT.cxx:478
 TVectorT.cxx:479
 TVectorT.cxx:480
 TVectorT.cxx:481
 TVectorT.cxx:482
 TVectorT.cxx:483
 TVectorT.cxx:484
 TVectorT.cxx:485
 TVectorT.cxx:486
 TVectorT.cxx:487
 TVectorT.cxx:488
 TVectorT.cxx:489
 TVectorT.cxx:490
 TVectorT.cxx:491
 TVectorT.cxx:492
 TVectorT.cxx:493
 TVectorT.cxx:494
 TVectorT.cxx:495
 TVectorT.cxx:496
 TVectorT.cxx:497
 TVectorT.cxx:498
 TVectorT.cxx:499
 TVectorT.cxx:500
 TVectorT.cxx:501
 TVectorT.cxx:502
 TVectorT.cxx:503
 TVectorT.cxx:504
 TVectorT.cxx:505
 TVectorT.cxx:506
 TVectorT.cxx:507
 TVectorT.cxx:508
 TVectorT.cxx:509
 TVectorT.cxx:510
 TVectorT.cxx:511
 TVectorT.cxx:512
 TVectorT.cxx:513
 TVectorT.cxx:514
 TVectorT.cxx:515
 TVectorT.cxx:516
 TVectorT.cxx:517
 TVectorT.cxx:518
 TVectorT.cxx:519
 TVectorT.cxx:520
 TVectorT.cxx:521
 TVectorT.cxx:522
 TVectorT.cxx:523
 TVectorT.cxx:524
 TVectorT.cxx:525
 TVectorT.cxx:526
 TVectorT.cxx:527
 TVectorT.cxx:528
 TVectorT.cxx:529
 TVectorT.cxx:530
 TVectorT.cxx:531
 TVectorT.cxx:532
 TVectorT.cxx:533
 TVectorT.cxx:534
 TVectorT.cxx:535
 TVectorT.cxx:536
 TVectorT.cxx:537
 TVectorT.cxx:538
 TVectorT.cxx:539
 TVectorT.cxx:540
 TVectorT.cxx:541
 TVectorT.cxx:542
 TVectorT.cxx:543
 TVectorT.cxx:544
 TVectorT.cxx:545
 TVectorT.cxx:546
 TVectorT.cxx:547
 TVectorT.cxx:548
 TVectorT.cxx:549
 TVectorT.cxx:550
 TVectorT.cxx:551
 TVectorT.cxx:552
 TVectorT.cxx:553
 TVectorT.cxx:554
 TVectorT.cxx:555
 TVectorT.cxx:556
 TVectorT.cxx:557
 TVectorT.cxx:558
 TVectorT.cxx:559
 TVectorT.cxx:560
 TVectorT.cxx:561
 TVectorT.cxx:562
 TVectorT.cxx:563
 TVectorT.cxx:564
 TVectorT.cxx:565
 TVectorT.cxx:566
 TVectorT.cxx:567
 TVectorT.cxx:568
 TVectorT.cxx:569
 TVectorT.cxx:570
 TVectorT.cxx:571
 TVectorT.cxx:572
 TVectorT.cxx:573
 TVectorT.cxx:574
 TVectorT.cxx:575
 TVectorT.cxx:576
 TVectorT.cxx:577
 TVectorT.cxx:578
 TVectorT.cxx:579
 TVectorT.cxx:580
 TVectorT.cxx:581
 TVectorT.cxx:582
 TVectorT.cxx:583
 TVectorT.cxx:584
 TVectorT.cxx:585
 TVectorT.cxx:586
 TVectorT.cxx:587
 TVectorT.cxx:588
 TVectorT.cxx:589
 TVectorT.cxx:590
 TVectorT.cxx:591
 TVectorT.cxx:592
 TVectorT.cxx:593
 TVectorT.cxx:594
 TVectorT.cxx:595
 TVectorT.cxx:596
 TVectorT.cxx:597
 TVectorT.cxx:598
 TVectorT.cxx:599
 TVectorT.cxx:600
 TVectorT.cxx:601
 TVectorT.cxx:602
 TVectorT.cxx:603
 TVectorT.cxx:604
 TVectorT.cxx:605
 TVectorT.cxx:606
 TVectorT.cxx:607
 TVectorT.cxx:608
 TVectorT.cxx:609
 TVectorT.cxx:610
 TVectorT.cxx:611
 TVectorT.cxx:612
 TVectorT.cxx:613
 TVectorT.cxx:614
 TVectorT.cxx:615
 TVectorT.cxx:616
 TVectorT.cxx:617
 TVectorT.cxx:618
 TVectorT.cxx:619
 TVectorT.cxx:620
 TVectorT.cxx:621
 TVectorT.cxx:622
 TVectorT.cxx:623
 TVectorT.cxx:624
 TVectorT.cxx:625
 TVectorT.cxx:626
 TVectorT.cxx:627
 TVectorT.cxx:628
 TVectorT.cxx:629
 TVectorT.cxx:630
 TVectorT.cxx:631
 TVectorT.cxx:632
 TVectorT.cxx:633
 TVectorT.cxx:634
 TVectorT.cxx:635
 TVectorT.cxx:636
 TVectorT.cxx:637
 TVectorT.cxx:638
 TVectorT.cxx:639
 TVectorT.cxx:640
 TVectorT.cxx:641
 TVectorT.cxx:642
 TVectorT.cxx:643
 TVectorT.cxx:644
 TVectorT.cxx:645
 TVectorT.cxx:646
 TVectorT.cxx:647
 TVectorT.cxx:648
 TVectorT.cxx:649
 TVectorT.cxx:650
 TVectorT.cxx:651
 TVectorT.cxx:652
 TVectorT.cxx:653
 TVectorT.cxx:654
 TVectorT.cxx:655
 TVectorT.cxx:656
 TVectorT.cxx:657
 TVectorT.cxx:658
 TVectorT.cxx:659
 TVectorT.cxx:660
 TVectorT.cxx:661
 TVectorT.cxx:662
 TVectorT.cxx:663
 TVectorT.cxx:664
 TVectorT.cxx:665
 TVectorT.cxx:666
 TVectorT.cxx:667
 TVectorT.cxx:668
 TVectorT.cxx:669
 TVectorT.cxx:670
 TVectorT.cxx:671
 TVectorT.cxx:672
 TVectorT.cxx:673
 TVectorT.cxx:674
 TVectorT.cxx:675
 TVectorT.cxx:676
 TVectorT.cxx:677
 TVectorT.cxx:678
 TVectorT.cxx:679
 TVectorT.cxx:680
 TVectorT.cxx:681
 TVectorT.cxx:682
 TVectorT.cxx:683
 TVectorT.cxx:684
 TVectorT.cxx:685
 TVectorT.cxx:686
 TVectorT.cxx:687
 TVectorT.cxx:688
 TVectorT.cxx:689
 TVectorT.cxx:690
 TVectorT.cxx:691
 TVectorT.cxx:692
 TVectorT.cxx:693
 TVectorT.cxx:694
 TVectorT.cxx:695
 TVectorT.cxx:696
 TVectorT.cxx:697
 TVectorT.cxx:698
 TVectorT.cxx:699
 TVectorT.cxx:700
 TVectorT.cxx:701
 TVectorT.cxx:702
 TVectorT.cxx:703
 TVectorT.cxx:704
 TVectorT.cxx:705
 TVectorT.cxx:706
 TVectorT.cxx:707
 TVectorT.cxx:708
 TVectorT.cxx:709
 TVectorT.cxx:710
 TVectorT.cxx:711
 TVectorT.cxx:712
 TVectorT.cxx:713
 TVectorT.cxx:714
 TVectorT.cxx:715
 TVectorT.cxx:716
 TVectorT.cxx:717
 TVectorT.cxx:718
 TVectorT.cxx:719
 TVectorT.cxx:720
 TVectorT.cxx:721
 TVectorT.cxx:722
 TVectorT.cxx:723
 TVectorT.cxx:724
 TVectorT.cxx:725
 TVectorT.cxx:726
 TVectorT.cxx:727
 TVectorT.cxx:728
 TVectorT.cxx:729
 TVectorT.cxx:730
 TVectorT.cxx:731
 TVectorT.cxx:732
 TVectorT.cxx:733
 TVectorT.cxx:734
 TVectorT.cxx:735
 TVectorT.cxx:736
 TVectorT.cxx:737
 TVectorT.cxx:738
 TVectorT.cxx:739
 TVectorT.cxx:740
 TVectorT.cxx:741
 TVectorT.cxx:742
 TVectorT.cxx:743
 TVectorT.cxx:744
 TVectorT.cxx:745
 TVectorT.cxx:746
 TVectorT.cxx:747
 TVectorT.cxx:748
 TVectorT.cxx:749
 TVectorT.cxx:750
 TVectorT.cxx:751
 TVectorT.cxx:752
 TVectorT.cxx:753
 TVectorT.cxx:754
 TVectorT.cxx:755
 TVectorT.cxx:756
 TVectorT.cxx:757
 TVectorT.cxx:758
 TVectorT.cxx:759
 TVectorT.cxx:760
 TVectorT.cxx:761
 TVectorT.cxx:762
 TVectorT.cxx:763
 TVectorT.cxx:764
 TVectorT.cxx:765
 TVectorT.cxx:766
 TVectorT.cxx:767
 TVectorT.cxx:768
 TVectorT.cxx:769
 TVectorT.cxx:770
 TVectorT.cxx:771
 TVectorT.cxx:772
 TVectorT.cxx:773
 TVectorT.cxx:774
 TVectorT.cxx:775
 TVectorT.cxx:776
 TVectorT.cxx:777
 TVectorT.cxx:778
 TVectorT.cxx:779
 TVectorT.cxx:780
 TVectorT.cxx:781
 TVectorT.cxx:782
 TVectorT.cxx:783
 TVectorT.cxx:784
 TVectorT.cxx:785
 TVectorT.cxx:786
 TVectorT.cxx:787
 TVectorT.cxx:788
 TVectorT.cxx:789
 TVectorT.cxx:790
 TVectorT.cxx:791
 TVectorT.cxx:792
 TVectorT.cxx:793
 TVectorT.cxx:794
 TVectorT.cxx:795
 TVectorT.cxx:796
 TVectorT.cxx:797
 TVectorT.cxx:798
 TVectorT.cxx:799
 TVectorT.cxx:800
 TVectorT.cxx:801
 TVectorT.cxx:802
 TVectorT.cxx:803
 TVectorT.cxx:804
 TVectorT.cxx:805
 TVectorT.cxx:806
 TVectorT.cxx:807
 TVectorT.cxx:808
 TVectorT.cxx:809
 TVectorT.cxx:810
 TVectorT.cxx:811
 TVectorT.cxx:812
 TVectorT.cxx:813
 TVectorT.cxx:814
 TVectorT.cxx:815
 TVectorT.cxx:816
 TVectorT.cxx:817
 TVectorT.cxx:818
 TVectorT.cxx:819
 TVectorT.cxx:820
 TVectorT.cxx:821
 TVectorT.cxx:822
 TVectorT.cxx:823
 TVectorT.cxx:824
 TVectorT.cxx:825
 TVectorT.cxx:826
 TVectorT.cxx:827
 TVectorT.cxx:828
 TVectorT.cxx:829
 TVectorT.cxx:830
 TVectorT.cxx:831
 TVectorT.cxx:832
 TVectorT.cxx:833
 TVectorT.cxx:834
 TVectorT.cxx:835
 TVectorT.cxx:836
 TVectorT.cxx:837
 TVectorT.cxx:838
 TVectorT.cxx:839
 TVectorT.cxx:840
 TVectorT.cxx:841
 TVectorT.cxx:842
 TVectorT.cxx:843
 TVectorT.cxx:844
 TVectorT.cxx:845
 TVectorT.cxx:846
 TVectorT.cxx:847
 TVectorT.cxx:848
 TVectorT.cxx:849
 TVectorT.cxx:850
 TVectorT.cxx:851
 TVectorT.cxx:852
 TVectorT.cxx:853
 TVectorT.cxx:854
 TVectorT.cxx:855
 TVectorT.cxx:856
 TVectorT.cxx:857
 TVectorT.cxx:858
 TVectorT.cxx:859
 TVectorT.cxx:860
 TVectorT.cxx:861
 TVectorT.cxx:862
 TVectorT.cxx:863
 TVectorT.cxx:864
 TVectorT.cxx:865
 TVectorT.cxx:866
 TVectorT.cxx:867
 TVectorT.cxx:868
 TVectorT.cxx:869
 TVectorT.cxx:870
 TVectorT.cxx:871
 TVectorT.cxx:872
 TVectorT.cxx:873
 TVectorT.cxx:874
 TVectorT.cxx:875
 TVectorT.cxx:876
 TVectorT.cxx:877
 TVectorT.cxx:878
 TVectorT.cxx:879
 TVectorT.cxx:880
 TVectorT.cxx:881
 TVectorT.cxx:882
 TVectorT.cxx:883
 TVectorT.cxx:884
 TVectorT.cxx:885
 TVectorT.cxx:886
 TVectorT.cxx:887
 TVectorT.cxx:888
 TVectorT.cxx:889
 TVectorT.cxx:890
 TVectorT.cxx:891
 TVectorT.cxx:892
 TVectorT.cxx:893
 TVectorT.cxx:894
 TVectorT.cxx:895
 TVectorT.cxx:896
 TVectorT.cxx:897
 TVectorT.cxx:898
 TVectorT.cxx:899
 TVectorT.cxx:900
 TVectorT.cxx:901
 TVectorT.cxx:902
 TVectorT.cxx:903
 TVectorT.cxx:904
 TVectorT.cxx:905
 TVectorT.cxx:906
 TVectorT.cxx:907
 TVectorT.cxx:908
 TVectorT.cxx:909
 TVectorT.cxx:910
 TVectorT.cxx:911
 TVectorT.cxx:912
 TVectorT.cxx:913
 TVectorT.cxx:914
 TVectorT.cxx:915
 TVectorT.cxx:916
 TVectorT.cxx:917
 TVectorT.cxx:918
 TVectorT.cxx:919
 TVectorT.cxx:920
 TVectorT.cxx:921
 TVectorT.cxx:922
 TVectorT.cxx:923
 TVectorT.cxx:924
 TVectorT.cxx:925
 TVectorT.cxx:926
 TVectorT.cxx:927
 TVectorT.cxx:928
 TVectorT.cxx:929
 TVectorT.cxx:930
 TVectorT.cxx:931
 TVectorT.cxx:932
 TVectorT.cxx:933
 TVectorT.cxx:934
 TVectorT.cxx:935
 TVectorT.cxx:936
 TVectorT.cxx:937
 TVectorT.cxx:938
 TVectorT.cxx:939
 TVectorT.cxx:940
 TVectorT.cxx:941
 TVectorT.cxx:942
 TVectorT.cxx:943
 TVectorT.cxx:944
 TVectorT.cxx:945
 TVectorT.cxx:946
 TVectorT.cxx:947
 TVectorT.cxx:948
 TVectorT.cxx:949
 TVectorT.cxx:950
 TVectorT.cxx:951
 TVectorT.cxx:952
 TVectorT.cxx:953
 TVectorT.cxx:954
 TVectorT.cxx:955
 TVectorT.cxx:956
 TVectorT.cxx:957
 TVectorT.cxx:958
 TVectorT.cxx:959
 TVectorT.cxx:960
 TVectorT.cxx:961
 TVectorT.cxx:962
 TVectorT.cxx:963
 TVectorT.cxx:964
 TVectorT.cxx:965
 TVectorT.cxx:966
 TVectorT.cxx:967
 TVectorT.cxx:968
 TVectorT.cxx:969
 TVectorT.cxx:970
 TVectorT.cxx:971
 TVectorT.cxx:972
 TVectorT.cxx:973
 TVectorT.cxx:974
 TVectorT.cxx:975
 TVectorT.cxx:976
 TVectorT.cxx:977
 TVectorT.cxx:978
 TVectorT.cxx:979
 TVectorT.cxx:980
 TVectorT.cxx:981
 TVectorT.cxx:982
 TVectorT.cxx:983
 TVectorT.cxx:984
 TVectorT.cxx:985
 TVectorT.cxx:986
 TVectorT.cxx:987
 TVectorT.cxx:988
 TVectorT.cxx:989
 TVectorT.cxx:990
 TVectorT.cxx:991
 TVectorT.cxx:992
 TVectorT.cxx:993
 TVectorT.cxx:994
 TVectorT.cxx:995
 TVectorT.cxx:996
 TVectorT.cxx:997
 TVectorT.cxx:998
 TVectorT.cxx:999
 TVectorT.cxx:1000
 TVectorT.cxx:1001
 TVectorT.cxx:1002
 TVectorT.cxx:1003
 TVectorT.cxx:1004
 TVectorT.cxx:1005
 TVectorT.cxx:1006
 TVectorT.cxx:1007
 TVectorT.cxx:1008
 TVectorT.cxx:1009
 TVectorT.cxx:1010
 TVectorT.cxx:1011
 TVectorT.cxx:1012
 TVectorT.cxx:1013
 TVectorT.cxx:1014
 TVectorT.cxx:1015
 TVectorT.cxx:1016
 TVectorT.cxx:1017
 TVectorT.cxx:1018
 TVectorT.cxx:1019
 TVectorT.cxx:1020
 TVectorT.cxx:1021
 TVectorT.cxx:1022
 TVectorT.cxx:1023
 TVectorT.cxx:1024
 TVectorT.cxx:1025
 TVectorT.cxx:1026
 TVectorT.cxx:1027
 TVectorT.cxx:1028
 TVectorT.cxx:1029
 TVectorT.cxx:1030
 TVectorT.cxx:1031
 TVectorT.cxx:1032
 TVectorT.cxx:1033
 TVectorT.cxx:1034
 TVectorT.cxx:1035
 TVectorT.cxx:1036
 TVectorT.cxx:1037
 TVectorT.cxx:1038
 TVectorT.cxx:1039
 TVectorT.cxx:1040
 TVectorT.cxx:1041
 TVectorT.cxx:1042
 TVectorT.cxx:1043
 TVectorT.cxx:1044
 TVectorT.cxx:1045
 TVectorT.cxx:1046
 TVectorT.cxx:1047
 TVectorT.cxx:1048
 TVectorT.cxx:1049
 TVectorT.cxx:1050
 TVectorT.cxx:1051
 TVectorT.cxx:1052
 TVectorT.cxx:1053
 TVectorT.cxx:1054
 TVectorT.cxx:1055
 TVectorT.cxx:1056
 TVectorT.cxx:1057
 TVectorT.cxx:1058
 TVectorT.cxx:1059
 TVectorT.cxx:1060
 TVectorT.cxx:1061
 TVectorT.cxx:1062
 TVectorT.cxx:1063
 TVectorT.cxx:1064
 TVectorT.cxx:1065
 TVectorT.cxx:1066
 TVectorT.cxx:1067
 TVectorT.cxx:1068
 TVectorT.cxx:1069
 TVectorT.cxx:1070
 TVectorT.cxx:1071
 TVectorT.cxx:1072
 TVectorT.cxx:1073
 TVectorT.cxx:1074
 TVectorT.cxx:1075
 TVectorT.cxx:1076
 TVectorT.cxx:1077
 TVectorT.cxx:1078
 TVectorT.cxx:1079
 TVectorT.cxx:1080
 TVectorT.cxx:1081
 TVectorT.cxx:1082
 TVectorT.cxx:1083
 TVectorT.cxx:1084
 TVectorT.cxx:1085
 TVectorT.cxx:1086
 TVectorT.cxx:1087
 TVectorT.cxx:1088
 TVectorT.cxx:1089
 TVectorT.cxx:1090
 TVectorT.cxx:1091
 TVectorT.cxx:1092
 TVectorT.cxx:1093
 TVectorT.cxx:1094
 TVectorT.cxx:1095
 TVectorT.cxx:1096
 TVectorT.cxx:1097
 TVectorT.cxx:1098
 TVectorT.cxx:1099
 TVectorT.cxx:1100
 TVectorT.cxx:1101
 TVectorT.cxx:1102
 TVectorT.cxx:1103
 TVectorT.cxx:1104
 TVectorT.cxx:1105
 TVectorT.cxx:1106
 TVectorT.cxx:1107
 TVectorT.cxx:1108
 TVectorT.cxx:1109
 TVectorT.cxx:1110
 TVectorT.cxx:1111
 TVectorT.cxx:1112
 TVectorT.cxx:1113
 TVectorT.cxx:1114
 TVectorT.cxx:1115
 TVectorT.cxx:1116
 TVectorT.cxx:1117
 TVectorT.cxx:1118
 TVectorT.cxx:1119
 TVectorT.cxx:1120
 TVectorT.cxx:1121
 TVectorT.cxx:1122
 TVectorT.cxx:1123
 TVectorT.cxx:1124
 TVectorT.cxx:1125
 TVectorT.cxx:1126
 TVectorT.cxx:1127
 TVectorT.cxx:1128
 TVectorT.cxx:1129
 TVectorT.cxx:1130
 TVectorT.cxx:1131
 TVectorT.cxx:1132
 TVectorT.cxx:1133
 TVectorT.cxx:1134
 TVectorT.cxx:1135
 TVectorT.cxx:1136
 TVectorT.cxx:1137
 TVectorT.cxx:1138
 TVectorT.cxx:1139
 TVectorT.cxx:1140
 TVectorT.cxx:1141
 TVectorT.cxx:1142
 TVectorT.cxx:1143
 TVectorT.cxx:1144
 TVectorT.cxx:1145
 TVectorT.cxx:1146
 TVectorT.cxx:1147
 TVectorT.cxx:1148
 TVectorT.cxx:1149
 TVectorT.cxx:1150
 TVectorT.cxx:1151
 TVectorT.cxx:1152
 TVectorT.cxx:1153
 TVectorT.cxx:1154
 TVectorT.cxx:1155
 TVectorT.cxx:1156
 TVectorT.cxx:1157
 TVectorT.cxx:1158
 TVectorT.cxx:1159
 TVectorT.cxx:1160
 TVectorT.cxx:1161
 TVectorT.cxx:1162
 TVectorT.cxx:1163
 TVectorT.cxx:1164
 TVectorT.cxx:1165
 TVectorT.cxx:1166
 TVectorT.cxx:1167
 TVectorT.cxx:1168
 TVectorT.cxx:1169
 TVectorT.cxx:1170
 TVectorT.cxx:1171
 TVectorT.cxx:1172
 TVectorT.cxx:1173
 TVectorT.cxx:1174
 TVectorT.cxx:1175
 TVectorT.cxx:1176
 TVectorT.cxx:1177
 TVectorT.cxx:1178
 TVectorT.cxx:1179
 TVectorT.cxx:1180
 TVectorT.cxx:1181
 TVectorT.cxx:1182
 TVectorT.cxx:1183
 TVectorT.cxx:1184
 TVectorT.cxx:1185
 TVectorT.cxx:1186
 TVectorT.cxx:1187
 TVectorT.cxx:1188
 TVectorT.cxx:1189
 TVectorT.cxx:1190
 TVectorT.cxx:1191
 TVectorT.cxx:1192
 TVectorT.cxx:1193
 TVectorT.cxx:1194
 TVectorT.cxx:1195
 TVectorT.cxx:1196
 TVectorT.cxx:1197
 TVectorT.cxx:1198
 TVectorT.cxx:1199
 TVectorT.cxx:1200
 TVectorT.cxx:1201
 TVectorT.cxx:1202
 TVectorT.cxx:1203
 TVectorT.cxx:1204
 TVectorT.cxx:1205
 TVectorT.cxx:1206
 TVectorT.cxx:1207
 TVectorT.cxx:1208
 TVectorT.cxx:1209
 TVectorT.cxx:1210
 TVectorT.cxx:1211
 TVectorT.cxx:1212
 TVectorT.cxx:1213
 TVectorT.cxx:1214
 TVectorT.cxx:1215
 TVectorT.cxx:1216
 TVectorT.cxx:1217
 TVectorT.cxx:1218
 TVectorT.cxx:1219
 TVectorT.cxx:1220
 TVectorT.cxx:1221
 TVectorT.cxx:1222
 TVectorT.cxx:1223
 TVectorT.cxx:1224
 TVectorT.cxx:1225
 TVectorT.cxx:1226
 TVectorT.cxx:1227
 TVectorT.cxx:1228
 TVectorT.cxx:1229
 TVectorT.cxx:1230
 TVectorT.cxx:1231
 TVectorT.cxx:1232
 TVectorT.cxx:1233
 TVectorT.cxx:1234
 TVectorT.cxx:1235
 TVectorT.cxx:1236
 TVectorT.cxx:1237
 TVectorT.cxx:1238
 TVectorT.cxx:1239
 TVectorT.cxx:1240
 TVectorT.cxx:1241
 TVectorT.cxx:1242
 TVectorT.cxx:1243
 TVectorT.cxx:1244
 TVectorT.cxx:1245
 TVectorT.cxx:1246
 TVectorT.cxx:1247
 TVectorT.cxx:1248
 TVectorT.cxx:1249
 TVectorT.cxx:1250
 TVectorT.cxx:1251
 TVectorT.cxx:1252
 TVectorT.cxx:1253
 TVectorT.cxx:1254
 TVectorT.cxx:1255
 TVectorT.cxx:1256
 TVectorT.cxx:1257
 TVectorT.cxx:1258
 TVectorT.cxx:1259
 TVectorT.cxx:1260
 TVectorT.cxx:1261
 TVectorT.cxx:1262
 TVectorT.cxx:1263
 TVectorT.cxx:1264
 TVectorT.cxx:1265
 TVectorT.cxx:1266
 TVectorT.cxx:1267
 TVectorT.cxx:1268
 TVectorT.cxx:1269
 TVectorT.cxx:1270
 TVectorT.cxx:1271
 TVectorT.cxx:1272
 TVectorT.cxx:1273
 TVectorT.cxx:1274
 TVectorT.cxx:1275
 TVectorT.cxx:1276
 TVectorT.cxx:1277
 TVectorT.cxx:1278
 TVectorT.cxx:1279
 TVectorT.cxx:1280
 TVectorT.cxx:1281
 TVectorT.cxx:1282
 TVectorT.cxx:1283
 TVectorT.cxx:1284
 TVectorT.cxx:1285
 TVectorT.cxx:1286
 TVectorT.cxx:1287
 TVectorT.cxx:1288
 TVectorT.cxx:1289
 TVectorT.cxx:1290
 TVectorT.cxx:1291
 TVectorT.cxx:1292
 TVectorT.cxx:1293
 TVectorT.cxx:1294
 TVectorT.cxx:1295
 TVectorT.cxx:1296
 TVectorT.cxx:1297
 TVectorT.cxx:1298
 TVectorT.cxx:1299
 TVectorT.cxx:1300
 TVectorT.cxx:1301
 TVectorT.cxx:1302
 TVectorT.cxx:1303
 TVectorT.cxx:1304
 TVectorT.cxx:1305
 TVectorT.cxx:1306
 TVectorT.cxx:1307
 TVectorT.cxx:1308
 TVectorT.cxx:1309
 TVectorT.cxx:1310
 TVectorT.cxx:1311
 TVectorT.cxx:1312
 TVectorT.cxx:1313
 TVectorT.cxx:1314
 TVectorT.cxx:1315
 TVectorT.cxx:1316
 TVectorT.cxx:1317
 TVectorT.cxx:1318
 TVectorT.cxx:1319
 TVectorT.cxx:1320
 TVectorT.cxx:1321
 TVectorT.cxx:1322
 TVectorT.cxx:1323
 TVectorT.cxx:1324
 TVectorT.cxx:1325
 TVectorT.cxx:1326
 TVectorT.cxx:1327
 TVectorT.cxx:1328
 TVectorT.cxx:1329
 TVectorT.cxx:1330
 TVectorT.cxx:1331
 TVectorT.cxx:1332
 TVectorT.cxx:1333
 TVectorT.cxx:1334
 TVectorT.cxx:1335
 TVectorT.cxx:1336
 TVectorT.cxx:1337
 TVectorT.cxx:1338
 TVectorT.cxx:1339
 TVectorT.cxx:1340
 TVectorT.cxx:1341
 TVectorT.cxx:1342
 TVectorT.cxx:1343
 TVectorT.cxx:1344
 TVectorT.cxx:1345
 TVectorT.cxx:1346
 TVectorT.cxx:1347
 TVectorT.cxx:1348
 TVectorT.cxx:1349
 TVectorT.cxx:1350
 TVectorT.cxx:1351
 TVectorT.cxx:1352
 TVectorT.cxx:1353
 TVectorT.cxx:1354
 TVectorT.cxx:1355
 TVectorT.cxx:1356
 TVectorT.cxx:1357
 TVectorT.cxx:1358
 TVectorT.cxx:1359
 TVectorT.cxx:1360
 TVectorT.cxx:1361
 TVectorT.cxx:1362
 TVectorT.cxx:1363
 TVectorT.cxx:1364
 TVectorT.cxx:1365
 TVectorT.cxx:1366
 TVectorT.cxx:1367
 TVectorT.cxx:1368
 TVectorT.cxx:1369
 TVectorT.cxx:1370
 TVectorT.cxx:1371
 TVectorT.cxx:1372
 TVectorT.cxx:1373
 TVectorT.cxx:1374
 TVectorT.cxx:1375
 TVectorT.cxx:1376
 TVectorT.cxx:1377
 TVectorT.cxx:1378
 TVectorT.cxx:1379
 TVectorT.cxx:1380
 TVectorT.cxx:1381
 TVectorT.cxx:1382
 TVectorT.cxx:1383
 TVectorT.cxx:1384
 TVectorT.cxx:1385
 TVectorT.cxx:1386
 TVectorT.cxx:1387
 TVectorT.cxx:1388
 TVectorT.cxx:1389
 TVectorT.cxx:1390
 TVectorT.cxx:1391
 TVectorT.cxx:1392
 TVectorT.cxx:1393
 TVectorT.cxx:1394
 TVectorT.cxx:1395
 TVectorT.cxx:1396
 TVectorT.cxx:1397
 TVectorT.cxx:1398
 TVectorT.cxx:1399
 TVectorT.cxx:1400
 TVectorT.cxx:1401
 TVectorT.cxx:1402
 TVectorT.cxx:1403
 TVectorT.cxx:1404
 TVectorT.cxx:1405
 TVectorT.cxx:1406
 TVectorT.cxx:1407
 TVectorT.cxx:1408
 TVectorT.cxx:1409
 TVectorT.cxx:1410
 TVectorT.cxx:1411
 TVectorT.cxx:1412
 TVectorT.cxx:1413
 TVectorT.cxx:1414
 TVectorT.cxx:1415
 TVectorT.cxx:1416
 TVectorT.cxx:1417
 TVectorT.cxx:1418
 TVectorT.cxx:1419
 TVectorT.cxx:1420
 TVectorT.cxx:1421
 TVectorT.cxx:1422
 TVectorT.cxx:1423
 TVectorT.cxx:1424
 TVectorT.cxx:1425
 TVectorT.cxx:1426
 TVectorT.cxx:1427
 TVectorT.cxx:1428
 TVectorT.cxx:1429
 TVectorT.cxx:1430
 TVectorT.cxx:1431
 TVectorT.cxx:1432
 TVectorT.cxx:1433
 TVectorT.cxx:1434
 TVectorT.cxx:1435
 TVectorT.cxx:1436
 TVectorT.cxx:1437
 TVectorT.cxx:1438
 TVectorT.cxx:1439
 TVectorT.cxx:1440
 TVectorT.cxx:1441
 TVectorT.cxx:1442
 TVectorT.cxx:1443
 TVectorT.cxx:1444
 TVectorT.cxx:1445
 TVectorT.cxx:1446
 TVectorT.cxx:1447
 TVectorT.cxx:1448
 TVectorT.cxx:1449
 TVectorT.cxx:1450
 TVectorT.cxx:1451
 TVectorT.cxx:1452
 TVectorT.cxx:1453
 TVectorT.cxx:1454
 TVectorT.cxx:1455
 TVectorT.cxx:1456
 TVectorT.cxx:1457
 TVectorT.cxx:1458
 TVectorT.cxx:1459
 TVectorT.cxx:1460
 TVectorT.cxx:1461
 TVectorT.cxx:1462
 TVectorT.cxx:1463
 TVectorT.cxx:1464
 TVectorT.cxx:1465
 TVectorT.cxx:1466
 TVectorT.cxx:1467
 TVectorT.cxx:1468
 TVectorT.cxx:1469
 TVectorT.cxx:1470
 TVectorT.cxx:1471
 TVectorT.cxx:1472
 TVectorT.cxx:1473
 TVectorT.cxx:1474
 TVectorT.cxx:1475
 TVectorT.cxx:1476
 TVectorT.cxx:1477
 TVectorT.cxx:1478
 TVectorT.cxx:1479
 TVectorT.cxx:1480
 TVectorT.cxx:1481
 TVectorT.cxx:1482
 TVectorT.cxx:1483
 TVectorT.cxx:1484
 TVectorT.cxx:1485
 TVectorT.cxx:1486
 TVectorT.cxx:1487
 TVectorT.cxx:1488
 TVectorT.cxx:1489
 TVectorT.cxx:1490
 TVectorT.cxx:1491
 TVectorT.cxx:1492
 TVectorT.cxx:1493
 TVectorT.cxx:1494
 TVectorT.cxx:1495
 TVectorT.cxx:1496
 TVectorT.cxx:1497
 TVectorT.cxx:1498
 TVectorT.cxx:1499
 TVectorT.cxx:1500
 TVectorT.cxx:1501
 TVectorT.cxx:1502
 TVectorT.cxx:1503
 TVectorT.cxx:1504
 TVectorT.cxx:1505
 TVectorT.cxx:1506
 TVectorT.cxx:1507
 TVectorT.cxx:1508
 TVectorT.cxx:1509
 TVectorT.cxx:1510
 TVectorT.cxx:1511
 TVectorT.cxx:1512
 TVectorT.cxx:1513
 TVectorT.cxx:1514
 TVectorT.cxx:1515
 TVectorT.cxx:1516
 TVectorT.cxx:1517
 TVectorT.cxx:1518
 TVectorT.cxx:1519
 TVectorT.cxx:1520
 TVectorT.cxx:1521
 TVectorT.cxx:1522
 TVectorT.cxx:1523
 TVectorT.cxx:1524
 TVectorT.cxx:1525
 TVectorT.cxx:1526
 TVectorT.cxx:1527
 TVectorT.cxx:1528
 TVectorT.cxx:1529
 TVectorT.cxx:1530
 TVectorT.cxx:1531
 TVectorT.cxx:1532
 TVectorT.cxx:1533
 TVectorT.cxx:1534
 TVectorT.cxx:1535
 TVectorT.cxx:1536
 TVectorT.cxx:1537
 TVectorT.cxx:1538
 TVectorT.cxx:1539
 TVectorT.cxx:1540
 TVectorT.cxx:1541
 TVectorT.cxx:1542
 TVectorT.cxx:1543
 TVectorT.cxx:1544
 TVectorT.cxx:1545
 TVectorT.cxx:1546
 TVectorT.cxx:1547
 TVectorT.cxx:1548
 TVectorT.cxx:1549
 TVectorT.cxx:1550
 TVectorT.cxx:1551
 TVectorT.cxx:1552
 TVectorT.cxx:1553
 TVectorT.cxx:1554
 TVectorT.cxx:1555
 TVectorT.cxx:1556
 TVectorT.cxx:1557
 TVectorT.cxx:1558
 TVectorT.cxx:1559
 TVectorT.cxx:1560
 TVectorT.cxx:1561
 TVectorT.cxx:1562
 TVectorT.cxx:1563
 TVectorT.cxx:1564
 TVectorT.cxx:1565
 TVectorT.cxx:1566
 TVectorT.cxx:1567
 TVectorT.cxx:1568
 TVectorT.cxx:1569
 TVectorT.cxx:1570
 TVectorT.cxx:1571
 TVectorT.cxx:1572
 TVectorT.cxx:1573
 TVectorT.cxx:1574
 TVectorT.cxx:1575
 TVectorT.cxx:1576
 TVectorT.cxx:1577
 TVectorT.cxx:1578
 TVectorT.cxx:1579
 TVectorT.cxx:1580
 TVectorT.cxx:1581
 TVectorT.cxx:1582
 TVectorT.cxx:1583
 TVectorT.cxx:1584
 TVectorT.cxx:1585
 TVectorT.cxx:1586
 TVectorT.cxx:1587
 TVectorT.cxx:1588
 TVectorT.cxx:1589
 TVectorT.cxx:1590
 TVectorT.cxx:1591
 TVectorT.cxx:1592
 TVectorT.cxx:1593
 TVectorT.cxx:1594
 TVectorT.cxx:1595
 TVectorT.cxx:1596
 TVectorT.cxx:1597
 TVectorT.cxx:1598
 TVectorT.cxx:1599
 TVectorT.cxx:1600
 TVectorT.cxx:1601
 TVectorT.cxx:1602
 TVectorT.cxx:1603
 TVectorT.cxx:1604
 TVectorT.cxx:1605
 TVectorT.cxx:1606
 TVectorT.cxx:1607
 TVectorT.cxx:1608
 TVectorT.cxx:1609
 TVectorT.cxx:1610
 TVectorT.cxx:1611
 TVectorT.cxx:1612
 TVectorT.cxx:1613
 TVectorT.cxx:1614
 TVectorT.cxx:1615
 TVectorT.cxx:1616
 TVectorT.cxx:1617
 TVectorT.cxx:1618
 TVectorT.cxx:1619
 TVectorT.cxx:1620
 TVectorT.cxx:1621
 TVectorT.cxx:1622
 TVectorT.cxx:1623
 TVectorT.cxx:1624
 TVectorT.cxx:1625
 TVectorT.cxx:1626
 TVectorT.cxx:1627
 TVectorT.cxx:1628
 TVectorT.cxx:1629
 TVectorT.cxx:1630
 TVectorT.cxx:1631
 TVectorT.cxx:1632
 TVectorT.cxx:1633
 TVectorT.cxx:1634
 TVectorT.cxx:1635
 TVectorT.cxx:1636
 TVectorT.cxx:1637
 TVectorT.cxx:1638
 TVectorT.cxx:1639
 TVectorT.cxx:1640
 TVectorT.cxx:1641
 TVectorT.cxx:1642
 TVectorT.cxx:1643
 TVectorT.cxx:1644
 TVectorT.cxx:1645
 TVectorT.cxx:1646
 TVectorT.cxx:1647
 TVectorT.cxx:1648
 TVectorT.cxx:1649
 TVectorT.cxx:1650
 TVectorT.cxx:1651
 TVectorT.cxx:1652
 TVectorT.cxx:1653
 TVectorT.cxx:1654
 TVectorT.cxx:1655
 TVectorT.cxx:1656
 TVectorT.cxx:1657
 TVectorT.cxx:1658
 TVectorT.cxx:1659
 TVectorT.cxx:1660
 TVectorT.cxx:1661
 TVectorT.cxx:1662
 TVectorT.cxx:1663
 TVectorT.cxx:1664
 TVectorT.cxx:1665
 TVectorT.cxx:1666
 TVectorT.cxx:1667
 TVectorT.cxx:1668
 TVectorT.cxx:1669
 TVectorT.cxx:1670
 TVectorT.cxx:1671
 TVectorT.cxx:1672
 TVectorT.cxx:1673
 TVectorT.cxx:1674
 TVectorT.cxx:1675
 TVectorT.cxx:1676
 TVectorT.cxx:1677
 TVectorT.cxx:1678
 TVectorT.cxx:1679
 TVectorT.cxx:1680
 TVectorT.cxx:1681
 TVectorT.cxx:1682
 TVectorT.cxx:1683
 TVectorT.cxx:1684
 TVectorT.cxx:1685
 TVectorT.cxx:1686
 TVectorT.cxx:1687
 TVectorT.cxx:1688
 TVectorT.cxx:1689
 TVectorT.cxx:1690
 TVectorT.cxx:1691
 TVectorT.cxx:1692
 TVectorT.cxx:1693
 TVectorT.cxx:1694
 TVectorT.cxx:1695
 TVectorT.cxx:1696
 TVectorT.cxx:1697
 TVectorT.cxx:1698
 TVectorT.cxx:1699
 TVectorT.cxx:1700
 TVectorT.cxx:1701
 TVectorT.cxx:1702
 TVectorT.cxx:1703
 TVectorT.cxx:1704
 TVectorT.cxx:1705
 TVectorT.cxx:1706
 TVectorT.cxx:1707
 TVectorT.cxx:1708
 TVectorT.cxx:1709
 TVectorT.cxx:1710
 TVectorT.cxx:1711
 TVectorT.cxx:1712
 TVectorT.cxx:1713
 TVectorT.cxx:1714
 TVectorT.cxx:1715
 TVectorT.cxx:1716
 TVectorT.cxx:1717
 TVectorT.cxx:1718
 TVectorT.cxx:1719
 TVectorT.cxx:1720
 TVectorT.cxx:1721
 TVectorT.cxx:1722
 TVectorT.cxx:1723
 TVectorT.cxx:1724
 TVectorT.cxx:1725
 TVectorT.cxx:1726
 TVectorT.cxx:1727
 TVectorT.cxx:1728
 TVectorT.cxx:1729
 TVectorT.cxx:1730
 TVectorT.cxx:1731
 TVectorT.cxx:1732
 TVectorT.cxx:1733
 TVectorT.cxx:1734
 TVectorT.cxx:1735
 TVectorT.cxx:1736
 TVectorT.cxx:1737
 TVectorT.cxx:1738
 TVectorT.cxx:1739
 TVectorT.cxx:1740
 TVectorT.cxx:1741
 TVectorT.cxx:1742
 TVectorT.cxx:1743
 TVectorT.cxx:1744
 TVectorT.cxx:1745
 TVectorT.cxx:1746
 TVectorT.cxx:1747
 TVectorT.cxx:1748
 TVectorT.cxx:1749
 TVectorT.cxx:1750
 TVectorT.cxx:1751
 TVectorT.cxx:1752
 TVectorT.cxx:1753
 TVectorT.cxx:1754
 TVectorT.cxx:1755
 TVectorT.cxx:1756
 TVectorT.cxx:1757
 TVectorT.cxx:1758
 TVectorT.cxx:1759
 TVectorT.cxx:1760
 TVectorT.cxx:1761
 TVectorT.cxx:1762
 TVectorT.cxx:1763
 TVectorT.cxx:1764
 TVectorT.cxx:1765
 TVectorT.cxx:1766
 TVectorT.cxx:1767
 TVectorT.cxx:1768
 TVectorT.cxx:1769
 TVectorT.cxx:1770
 TVectorT.cxx:1771
 TVectorT.cxx:1772
 TVectorT.cxx:1773
 TVectorT.cxx:1774
 TVectorT.cxx:1775
 TVectorT.cxx:1776
 TVectorT.cxx:1777
 TVectorT.cxx:1778
 TVectorT.cxx:1779
 TVectorT.cxx:1780
 TVectorT.cxx:1781
 TVectorT.cxx:1782
 TVectorT.cxx:1783
 TVectorT.cxx:1784
 TVectorT.cxx:1785
 TVectorT.cxx:1786
 TVectorT.cxx:1787
 TVectorT.cxx:1788
 TVectorT.cxx:1789
 TVectorT.cxx:1790
 TVectorT.cxx:1791
 TVectorT.cxx:1792
 TVectorT.cxx:1793
 TVectorT.cxx:1794
 TVectorT.cxx:1795
 TVectorT.cxx:1796
 TVectorT.cxx:1797
 TVectorT.cxx:1798
 TVectorT.cxx:1799
 TVectorT.cxx:1800
 TVectorT.cxx:1801
 TVectorT.cxx:1802
 TVectorT.cxx:1803
 TVectorT.cxx:1804
 TVectorT.cxx:1805
 TVectorT.cxx:1806
 TVectorT.cxx:1807
 TVectorT.cxx:1808
 TVectorT.cxx:1809
 TVectorT.cxx:1810
 TVectorT.cxx:1811
 TVectorT.cxx:1812
 TVectorT.cxx:1813
 TVectorT.cxx:1814
 TVectorT.cxx:1815
 TVectorT.cxx:1816
 TVectorT.cxx:1817
 TVectorT.cxx:1818
 TVectorT.cxx:1819
 TVectorT.cxx:1820
 TVectorT.cxx:1821
 TVectorT.cxx:1822
 TVectorT.cxx:1823
 TVectorT.cxx:1824
 TVectorT.cxx:1825
 TVectorT.cxx:1826
 TVectorT.cxx:1827
 TVectorT.cxx:1828
 TVectorT.cxx:1829
 TVectorT.cxx:1830
 TVectorT.cxx:1831
 TVectorT.cxx:1832
 TVectorT.cxx:1833
 TVectorT.cxx:1834
 TVectorT.cxx:1835
 TVectorT.cxx:1836
 TVectorT.cxx:1837
 TVectorT.cxx:1838
 TVectorT.cxx:1839
 TVectorT.cxx:1840
 TVectorT.cxx:1841
 TVectorT.cxx:1842
 TVectorT.cxx:1843
 TVectorT.cxx:1844
 TVectorT.cxx:1845
 TVectorT.cxx:1846
 TVectorT.cxx:1847
 TVectorT.cxx:1848
 TVectorT.cxx:1849
 TVectorT.cxx:1850
 TVectorT.cxx:1851
 TVectorT.cxx:1852
 TVectorT.cxx:1853
 TVectorT.cxx:1854
 TVectorT.cxx:1855
 TVectorT.cxx:1856
 TVectorT.cxx:1857
 TVectorT.cxx:1858
 TVectorT.cxx:1859
 TVectorT.cxx:1860
 TVectorT.cxx:1861
 TVectorT.cxx:1862
 TVectorT.cxx:1863
 TVectorT.cxx:1864
 TVectorT.cxx:1865
 TVectorT.cxx:1866
 TVectorT.cxx:1867
 TVectorT.cxx:1868
 TVectorT.cxx:1869
 TVectorT.cxx:1870
 TVectorT.cxx:1871
 TVectorT.cxx:1872
 TVectorT.cxx:1873
 TVectorT.cxx:1874
 TVectorT.cxx:1875
 TVectorT.cxx:1876
 TVectorT.cxx:1877
 TVectorT.cxx:1878
 TVectorT.cxx:1879
 TVectorT.cxx:1880
 TVectorT.cxx:1881
 TVectorT.cxx:1882
 TVectorT.cxx:1883
 TVectorT.cxx:1884
 TVectorT.cxx:1885
 TVectorT.cxx:1886
 TVectorT.cxx:1887
 TVectorT.cxx:1888
 TVectorT.cxx:1889
 TVectorT.cxx:1890
 TVectorT.cxx:1891
 TVectorT.cxx:1892
 TVectorT.cxx:1893
 TVectorT.cxx:1894
 TVectorT.cxx:1895
 TVectorT.cxx:1896
 TVectorT.cxx:1897
 TVectorT.cxx:1898
 TVectorT.cxx:1899
 TVectorT.cxx:1900
 TVectorT.cxx:1901
 TVectorT.cxx:1902
 TVectorT.cxx:1903
 TVectorT.cxx:1904
 TVectorT.cxx:1905
 TVectorT.cxx:1906
 TVectorT.cxx:1907
 TVectorT.cxx:1908
 TVectorT.cxx:1909
 TVectorT.cxx:1910
 TVectorT.cxx:1911
 TVectorT.cxx:1912
 TVectorT.cxx:1913
 TVectorT.cxx:1914
 TVectorT.cxx:1915
 TVectorT.cxx:1916
 TVectorT.cxx:1917
 TVectorT.cxx:1918
 TVectorT.cxx:1919
 TVectorT.cxx:1920
 TVectorT.cxx:1921
 TVectorT.cxx:1922
 TVectorT.cxx:1923
 TVectorT.cxx:1924
 TVectorT.cxx:1925
 TVectorT.cxx:1926
 TVectorT.cxx:1927
 TVectorT.cxx:1928
 TVectorT.cxx:1929
 TVectorT.cxx:1930
 TVectorT.cxx:1931
 TVectorT.cxx:1932
 TVectorT.cxx:1933
 TVectorT.cxx:1934
 TVectorT.cxx:1935
 TVectorT.cxx:1936
 TVectorT.cxx:1937
 TVectorT.cxx:1938
 TVectorT.cxx:1939
 TVectorT.cxx:1940
 TVectorT.cxx:1941
 TVectorT.cxx:1942
 TVectorT.cxx:1943
 TVectorT.cxx:1944
 TVectorT.cxx:1945
 TVectorT.cxx:1946
 TVectorT.cxx:1947
 TVectorT.cxx:1948
 TVectorT.cxx:1949
 TVectorT.cxx:1950
 TVectorT.cxx:1951
 TVectorT.cxx:1952
 TVectorT.cxx:1953
 TVectorT.cxx:1954
 TVectorT.cxx:1955
 TVectorT.cxx:1956
 TVectorT.cxx:1957
 TVectorT.cxx:1958
 TVectorT.cxx:1959
 TVectorT.cxx:1960
 TVectorT.cxx:1961
 TVectorT.cxx:1962
 TVectorT.cxx:1963
 TVectorT.cxx:1964
 TVectorT.cxx:1965
 TVectorT.cxx:1966
 TVectorT.cxx:1967
 TVectorT.cxx:1968
 TVectorT.cxx:1969
 TVectorT.cxx:1970
 TVectorT.cxx:1971
 TVectorT.cxx:1972
 TVectorT.cxx:1973
 TVectorT.cxx:1974
 TVectorT.cxx:1975
 TVectorT.cxx:1976
 TVectorT.cxx:1977
 TVectorT.cxx:1978
 TVectorT.cxx:1979
 TVectorT.cxx:1980
 TVectorT.cxx:1981
 TVectorT.cxx:1982
 TVectorT.cxx:1983
 TVectorT.cxx:1984
 TVectorT.cxx:1985
 TVectorT.cxx:1986
 TVectorT.cxx:1987
 TVectorT.cxx:1988
 TVectorT.cxx:1989
 TVectorT.cxx:1990
 TVectorT.cxx:1991
 TVectorT.cxx:1992
 TVectorT.cxx:1993
 TVectorT.cxx:1994
 TVectorT.cxx:1995
 TVectorT.cxx:1996
 TVectorT.cxx:1997
 TVectorT.cxx:1998
 TVectorT.cxx:1999
 TVectorT.cxx:2000
 TVectorT.cxx:2001
 TVectorT.cxx:2002
 TVectorT.cxx:2003
 TVectorT.cxx:2004
 TVectorT.cxx:2005
 TVectorT.cxx:2006
 TVectorT.cxx:2007
 TVectorT.cxx:2008
 TVectorT.cxx:2009
 TVectorT.cxx:2010
 TVectorT.cxx:2011
 TVectorT.cxx:2012
 TVectorT.cxx:2013
 TVectorT.cxx:2014
 TVectorT.cxx:2015
 TVectorT.cxx:2016
 TVectorT.cxx:2017
 TVectorT.cxx:2018
 TVectorT.cxx:2019
 TVectorT.cxx:2020
 TVectorT.cxx:2021
 TVectorT.cxx:2022
 TVectorT.cxx:2023
 TVectorT.cxx:2024
 TVectorT.cxx:2025
 TVectorT.cxx:2026
 TVectorT.cxx:2027
 TVectorT.cxx:2028
 TVectorT.cxx:2029
 TVectorT.cxx:2030
 TVectorT.cxx:2031
 TVectorT.cxx:2032
 TVectorT.cxx:2033
 TVectorT.cxx:2034
 TVectorT.cxx:2035
 TVectorT.cxx:2036
 TVectorT.cxx:2037
 TVectorT.cxx:2038
 TVectorT.cxx:2039
 TVectorT.cxx:2040
 TVectorT.cxx:2041
 TVectorT.cxx:2042
 TVectorT.cxx:2043
 TVectorT.cxx:2044
 TVectorT.cxx:2045
 TVectorT.cxx:2046
 TVectorT.cxx:2047
 TVectorT.cxx:2048
 TVectorT.cxx:2049
 TVectorT.cxx:2050
 TVectorT.cxx:2051
 TVectorT.cxx:2052
 TVectorT.cxx:2053
 TVectorT.cxx:2054
 TVectorT.cxx:2055
 TVectorT.cxx:2056
 TVectorT.cxx:2057
 TVectorT.cxx:2058
 TVectorT.cxx:2059
 TVectorT.cxx:2060
 TVectorT.cxx:2061
 TVectorT.cxx:2062
 TVectorT.cxx:2063
 TVectorT.cxx:2064
 TVectorT.cxx:2065
 TVectorT.cxx:2066
 TVectorT.cxx:2067
 TVectorT.cxx:2068
 TVectorT.cxx:2069
 TVectorT.cxx:2070
 TVectorT.cxx:2071
 TVectorT.cxx:2072
 TVectorT.cxx:2073
 TVectorT.cxx:2074
 TVectorT.cxx:2075
 TVectorT.cxx:2076
 TVectorT.cxx:2077
 TVectorT.cxx:2078
 TVectorT.cxx:2079
 TVectorT.cxx:2080
 TVectorT.cxx:2081
 TVectorT.cxx:2082
 TVectorT.cxx:2083
 TVectorT.cxx:2084
 TVectorT.cxx:2085
 TVectorT.cxx:2086
 TVectorT.cxx:2087
 TVectorT.cxx:2088
 TVectorT.cxx:2089
 TVectorT.cxx:2090
 TVectorT.cxx:2091
 TVectorT.cxx:2092
 TVectorT.cxx:2093
 TVectorT.cxx:2094
 TVectorT.cxx:2095
 TVectorT.cxx:2096
 TVectorT.cxx:2097
 TVectorT.cxx:2098
 TVectorT.cxx:2099
 TVectorT.cxx:2100
 TVectorT.cxx:2101
 TVectorT.cxx:2102
 TVectorT.cxx:2103
 TVectorT.cxx:2104
 TVectorT.cxx:2105
 TVectorT.cxx:2106
 TVectorT.cxx:2107
 TVectorT.cxx:2108
 TVectorT.cxx:2109
 TVectorT.cxx:2110
 TVectorT.cxx:2111
 TVectorT.cxx:2112
 TVectorT.cxx:2113
 TVectorT.cxx:2114
 TVectorT.cxx:2115
 TVectorT.cxx:2116
 TVectorT.cxx:2117
 TVectorT.cxx:2118
 TVectorT.cxx:2119
 TVectorT.cxx:2120
 TVectorT.cxx:2121
 TVectorT.cxx:2122
 TVectorT.cxx:2123
 TVectorT.cxx:2124
 TVectorT.cxx:2125
 TVectorT.cxx:2126
 TVectorT.cxx:2127
 TVectorT.cxx:2128
 TVectorT.cxx:2129
 TVectorT.cxx:2130
 TVectorT.cxx:2131
 TVectorT.cxx:2132
 TVectorT.cxx:2133
 TVectorT.cxx:2134
 TVectorT.cxx:2135
 TVectorT.cxx:2136
 TVectorT.cxx:2137
 TVectorT.cxx:2138
 TVectorT.cxx:2139
 TVectorT.cxx:2140
 TVectorT.cxx:2141
 TVectorT.cxx:2142
 TVectorT.cxx:2143
 TVectorT.cxx:2144
 TVectorT.cxx:2145
 TVectorT.cxx:2146
 TVectorT.cxx:2147
 TVectorT.cxx:2148
 TVectorT.cxx:2149
 TVectorT.cxx:2150
 TVectorT.cxx:2151
 TVectorT.cxx:2152
 TVectorT.cxx:2153
 TVectorT.cxx:2154
 TVectorT.cxx:2155
 TVectorT.cxx:2156
 TVectorT.cxx:2157
 TVectorT.cxx:2158
 TVectorT.cxx:2159
 TVectorT.cxx:2160
 TVectorT.cxx:2161
 TVectorT.cxx:2162
 TVectorT.cxx:2163
 TVectorT.cxx:2164
 TVectorT.cxx:2165
 TVectorT.cxx:2166
 TVectorT.cxx:2167
 TVectorT.cxx:2168
 TVectorT.cxx:2169
 TVectorT.cxx:2170
 TVectorT.cxx:2171
 TVectorT.cxx:2172
 TVectorT.cxx:2173
 TVectorT.cxx:2174
 TVectorT.cxx:2175
 TVectorT.cxx:2176
 TVectorT.cxx:2177
 TVectorT.cxx:2178
 TVectorT.cxx:2179
 TVectorT.cxx:2180
 TVectorT.cxx:2181
 TVectorT.cxx:2182
 TVectorT.cxx:2183
 TVectorT.cxx:2184
 TVectorT.cxx:2185
 TVectorT.cxx:2186
 TVectorT.cxx:2187
 TVectorT.cxx:2188
 TVectorT.cxx:2189
 TVectorT.cxx:2190
 TVectorT.cxx:2191
 TVectorT.cxx:2192
 TVectorT.cxx:2193
 TVectorT.cxx:2194
 TVectorT.cxx:2195
 TVectorT.cxx:2196
 TVectorT.cxx:2197
 TVectorT.cxx:2198
 TVectorT.cxx:2199
 TVectorT.cxx:2200
 TVectorT.cxx:2201
 TVectorT.cxx:2202
 TVectorT.cxx:2203
 TVectorT.cxx:2204
 TVectorT.cxx:2205
 TVectorT.cxx:2206
 TVectorT.cxx:2207
 TVectorT.cxx:2208
 TVectorT.cxx:2209
 TVectorT.cxx:2210
 TVectorT.cxx:2211
 TVectorT.cxx:2212
 TVectorT.cxx:2213
 TVectorT.cxx:2214
 TVectorT.cxx:2215
 TVectorT.cxx:2216
 TVectorT.cxx:2217
 TVectorT.cxx:2218
 TVectorT.cxx:2219
 TVectorT.cxx:2220
 TVectorT.cxx:2221
 TVectorT.cxx:2222
 TVectorT.cxx:2223
 TVectorT.cxx:2224
 TVectorT.cxx:2225
 TVectorT.cxx:2226
 TVectorT.cxx:2227
 TVectorT.cxx:2228
 TVectorT.cxx:2229
 TVectorT.cxx:2230
 TVectorT.cxx:2231
 TVectorT.cxx:2232
 TVectorT.cxx:2233
 TVectorT.cxx:2234
 TVectorT.cxx:2235
 TVectorT.cxx:2236
 TVectorT.cxx:2237
 TVectorT.cxx:2238
 TVectorT.cxx:2239
 TVectorT.cxx:2240
 TVectorT.cxx:2241
 TVectorT.cxx:2242
 TVectorT.cxx:2243
 TVectorT.cxx:2244
 TVectorT.cxx:2245
 TVectorT.cxx:2246
 TVectorT.cxx:2247
 TVectorT.cxx:2248
 TVectorT.cxx:2249
 TVectorT.cxx:2250
 TVectorT.cxx:2251
 TVectorT.cxx:2252
 TVectorT.cxx:2253
 TVectorT.cxx:2254
 TVectorT.cxx:2255
 TVectorT.cxx:2256
 TVectorT.cxx:2257
 TVectorT.cxx:2258
 TVectorT.cxx:2259
 TVectorT.cxx:2260
 TVectorT.cxx:2261
 TVectorT.cxx:2262
 TVectorT.cxx:2263
 TVectorT.cxx:2264
 TVectorT.cxx:2265
 TVectorT.cxx:2266
 TVectorT.cxx:2267
 TVectorT.cxx:2268
 TVectorT.cxx:2269
 TVectorT.cxx:2270
 TVectorT.cxx:2271
 TVectorT.cxx:2272
 TVectorT.cxx:2273
 TVectorT.cxx:2274
 TVectorT.cxx:2275
 TVectorT.cxx:2276
 TVectorT.cxx:2277
 TVectorT.cxx:2278
 TVectorT.cxx:2279
 TVectorT.cxx:2280
 TVectorT.cxx:2281
 TVectorT.cxx:2282
 TVectorT.cxx:2283
 TVectorT.cxx:2284
 TVectorT.cxx:2285
 TVectorT.cxx:2286
 TVectorT.cxx:2287
 TVectorT.cxx:2288
 TVectorT.cxx:2289
 TVectorT.cxx:2290
 TVectorT.cxx:2291
 TVectorT.cxx:2292
 TVectorT.cxx:2293
 TVectorT.cxx:2294
 TVectorT.cxx:2295
 TVectorT.cxx:2296
 TVectorT.cxx:2297
 TVectorT.cxx:2298
 TVectorT.cxx:2299
 TVectorT.cxx:2300
 TVectorT.cxx:2301
 TVectorT.cxx:2302
 TVectorT.cxx:2303
 TVectorT.cxx:2304
 TVectorT.cxx:2305
 TVectorT.cxx:2306
 TVectorT.cxx:2307
 TVectorT.cxx:2308
 TVectorT.cxx:2309
 TVectorT.cxx:2310
 TVectorT.cxx:2311
 TVectorT.cxx:2312
 TVectorT.cxx:2313
 TVectorT.cxx:2314
 TVectorT.cxx:2315
 TVectorT.cxx:2316
 TVectorT.cxx:2317
 TVectorT.cxx:2318
 TVectorT.cxx:2319
 TVectorT.cxx:2320
 TVectorT.cxx:2321
 TVectorT.cxx:2322
 TVectorT.cxx:2323
 TVectorT.cxx:2324
 TVectorT.cxx:2325
 TVectorT.cxx:2326
 TVectorT.cxx:2327
 TVectorT.cxx:2328
 TVectorT.cxx:2329
 TVectorT.cxx:2330
 TVectorT.cxx:2331
 TVectorT.cxx:2332
 TVectorT.cxx:2333
 TVectorT.cxx:2334
 TVectorT.cxx:2335
 TVectorT.cxx:2336
 TVectorT.cxx:2337
 TVectorT.cxx:2338
 TVectorT.cxx:2339
 TVectorT.cxx:2340
 TVectorT.cxx:2341
 TVectorT.cxx:2342
 TVectorT.cxx:2343
 TVectorT.cxx:2344
 TVectorT.cxx:2345
 TVectorT.cxx:2346
 TVectorT.cxx:2347
 TVectorT.cxx:2348
 TVectorT.cxx:2349
 TVectorT.cxx:2350
 TVectorT.cxx:2351
 TVectorT.cxx:2352
 TVectorT.cxx:2353
 TVectorT.cxx:2354
 TVectorT.cxx:2355
 TVectorT.cxx:2356
 TVectorT.cxx:2357
 TVectorT.cxx:2358
 TVectorT.cxx:2359
 TVectorT.cxx:2360
 TVectorT.cxx:2361
 TVectorT.cxx:2362
 TVectorT.cxx:2363
 TVectorT.cxx:2364
 TVectorT.cxx:2365
 TVectorT.cxx:2366
 TVectorT.cxx:2367
 TVectorT.cxx:2368
 TVectorT.cxx:2369
 TVectorT.cxx:2370
 TVectorT.cxx:2371
 TVectorT.cxx:2372
 TVectorT.cxx:2373
 TVectorT.cxx:2374
 TVectorT.cxx:2375
 TVectorT.cxx:2376
 TVectorT.cxx:2377
 TVectorT.cxx:2378
 TVectorT.cxx:2379
 TVectorT.cxx:2380
 TVectorT.cxx:2381
 TVectorT.cxx:2382
 TVectorT.cxx:2383
 TVectorT.cxx:2384
 TVectorT.cxx:2385
 TVectorT.cxx:2386
 TVectorT.cxx:2387
 TVectorT.cxx:2388
 TVectorT.cxx:2389
 TVectorT.cxx:2390
 TVectorT.cxx:2391
 TVectorT.cxx:2392
 TVectorT.cxx:2393
 TVectorT.cxx:2394
 TVectorT.cxx:2395
 TVectorT.cxx:2396
 TVectorT.cxx:2397
 TVectorT.cxx:2398
 TVectorT.cxx:2399
 TVectorT.cxx:2400
 TVectorT.cxx:2401
 TVectorT.cxx:2402
 TVectorT.cxx:2403
 TVectorT.cxx:2404
 TVectorT.cxx:2405
 TVectorT.cxx:2406
 TVectorT.cxx:2407
 TVectorT.cxx:2408
 TVectorT.cxx:2409
 TVectorT.cxx:2410
 TVectorT.cxx:2411
 TVectorT.cxx:2412
 TVectorT.cxx:2413
 TVectorT.cxx:2414
 TVectorT.cxx:2415
 TVectorT.cxx:2416
 TVectorT.cxx:2417
 TVectorT.cxx:2418
 TVectorT.cxx:2419
 TVectorT.cxx:2420
 TVectorT.cxx:2421
 TVectorT.cxx:2422
 TVectorT.cxx:2423
 TVectorT.cxx:2424
 TVectorT.cxx:2425
 TVectorT.cxx:2426
 TVectorT.cxx:2427
 TVectorT.cxx:2428
 TVectorT.cxx:2429
 TVectorT.cxx:2430
 TVectorT.cxx:2431
 TVectorT.cxx:2432
 TVectorT.cxx:2433
 TVectorT.cxx:2434
 TVectorT.cxx:2435
 TVectorT.cxx:2436
 TVectorT.cxx:2437
 TVectorT.cxx:2438
 TVectorT.cxx:2439
 TVectorT.cxx:2440
 TVectorT.cxx:2441
 TVectorT.cxx:2442
 TVectorT.cxx:2443
 TVectorT.cxx:2444
 TVectorT.cxx:2445
 TVectorT.cxx:2446
 TVectorT.cxx:2447
 TVectorT.cxx:2448
 TVectorT.cxx:2449
 TVectorT.cxx:2450
 TVectorT.cxx:2451
 TVectorT.cxx:2452
 TVectorT.cxx:2453
 TVectorT.cxx:2454
 TVectorT.cxx:2455
 TVectorT.cxx:2456
 TVectorT.cxx:2457
 TVectorT.cxx:2458
 TVectorT.cxx:2459
 TVectorT.cxx:2460
 TVectorT.cxx:2461
 TVectorT.cxx:2462
 TVectorT.cxx:2463
 TVectorT.cxx:2464
 TVectorT.cxx:2465
 TVectorT.cxx:2466
 TVectorT.cxx:2467
 TVectorT.cxx:2468
 TVectorT.cxx:2469
 TVectorT.cxx:2470
 TVectorT.cxx:2471
 TVectorT.cxx:2472
 TVectorT.cxx:2473
 TVectorT.cxx:2474
 TVectorT.cxx:2475