#include "TMatrixTBase.h"
#include "TVectorT.h"
#include "TClass.h"
Int_t gMatrixCheck = 1;
#ifndef R__ALPHA
templateClassImp(TMatrixTBase)
#endif
template<class Element>
void TMatrixTBase<Element>::DoubleLexSort(Int_t n,Int_t *first,Int_t *second,Element *data)
{
const int incs[] = {1,5,19,41,109,209,505,929,2161,3905,8929,16001,INT_MAX};
Int_t kinc = 0;
while (incs[kinc] <= n/2)
kinc++;
kinc -= 1;
for( ; kinc >= 0; kinc--) {
const Int_t inc = incs[kinc];
for (Int_t k = inc; k < n; k++) {
const Element tmp = data[k];
const Int_t fi = first [k];
const Int_t se = second[k];
Int_t j;
for (j = k; j >= inc; j -= inc) {
if ( fi < first[j-inc] || (fi == first[j-inc] && se < second[j-inc]) ) {
data [j] = data [j-inc];
first [j] = first [j-inc];
second[j] = second[j-inc];
} else
break;
}
data [j] = tmp;
first [j] = fi;
second[j] = se;
}
}
}
template<class Element>
void TMatrixTBase<Element>::IndexedLexSort(Int_t n,Int_t *first,Int_t swapFirst,
Int_t *second,Int_t swapSecond,Int_t *index)
{
const int incs[] = {1,5,19,41,109,209,505,929,2161,3905,8929,16001,INT_MAX};
Int_t kinc = 0;
while (incs[kinc] <= n/2)
kinc++;
kinc -= 1;
for( ; kinc >= 0; kinc--) {
const Int_t inc = incs[kinc];
if ( !swapFirst && !swapSecond ) {
for (Int_t k = inc; k < n; k++) {
const Int_t ktemp = index[k];
const Int_t fi = first [ktemp];
const Int_t se = second[ktemp];
Int_t j;
for (j = k; j >= inc; j -= inc) {
if (fi < first[index[j-inc]] || (fi == first[index[j-inc]] && se < second[index[j-inc]])) {
index[j] = index[j-inc];
} else {
break;
}
}
index[j] = ktemp;
}
} else if ( swapSecond && !swapFirst ) {
for (Int_t k = inc; k < n; k++) {
const Int_t ktemp = index[k];
const Int_t fi = first [ktemp];
const Int_t se = second[k];
Int_t j;
for (j = k; j >= inc; j -= inc) {
if (fi < first[index[j-inc]] || (fi == first[index[j-inc]] && se < second[j-inc])) {
index [j] = index[j-inc];
second[j] = second[j-inc];
} else {
break;
}
}
index[j] = ktemp;
second[j] = se;
}
} else if (swapFirst && !swapSecond) {
for (Int_t k = inc; k < n; k++ ) {
const Int_t ktemp = index[k];
const Int_t fi = first[k];
const Int_t se = second[ktemp];
Int_t j;
for (j = k; j >= inc; j -= inc) {
if ( fi < first[j-inc] || (fi == first[j-inc] && se < second[ index[j-inc]])) {
index[j] = index[j-inc];
first[j] = first[j-inc];
} else {
break;
}
}
index[j] = ktemp;
first[j] = fi;
}
} else {
for (Int_t k = inc; k < n; k++ ) {
const Int_t ktemp = index[k];
const Int_t fi = first [k];
const Int_t se = second[k];
Int_t j;
for (j = k; j >= inc; j -= inc) {
if ( fi < first[j-inc] || (fi == first[j-inc] && se < second[j-inc])) {
index [j] = index [j-inc];
first [j] = first [j-inc];
second[j] = second[j-inc];
} else {
break;
}
}
index[j] = ktemp;
first[j] = fi;
second[j] = se;
}
}
}
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::SetMatrixArray(const Element *data,Option_t *option)
{
R__ASSERT(IsValid());
TString opt = option;
opt.ToUpper();
Element *elem = GetMatrixArray();
if (opt.Contains("F")) {
for (Int_t irow = 0; irow < fNrows; irow++) {
const Int_t off1 = irow*fNcols;
Int_t off2 = 0;
for (Int_t icol = 0; icol < fNcols; icol++) {
elem[off1+icol] = data[off2+irow];
off2 += fNrows;
}
}
}
else
memcpy(elem,data,fNelems*sizeof(Element));
return *this;
}
template<class Element>
Bool_t TMatrixTBase<Element>::IsSymmetric() const
{
R__ASSERT(IsValid());
if ((fNrows != fNcols) || (fRowLwb != fColLwb))
return kFALSE;
const Element * const elem = GetMatrixArray();
for (Int_t irow = 0; irow < fNrows; irow++) {
const Int_t rowOff = irow*fNcols;
Int_t colOff = 0;
for (Int_t icol = 0; icol < irow; icol++) {
if (elem[rowOff+icol] != elem[colOff+irow])
return kFALSE;
colOff += fNrows;
}
}
return kTRUE;
}
template<class Element>
void TMatrixTBase<Element>::GetMatrix2Array(Element *data,Option_t *option) const
{
R__ASSERT(IsValid());
TString opt = option;
opt.ToUpper();
const Element * const elem = GetMatrixArray();
if (opt.Contains("F")) {
for (Int_t irow = 0; irow < fNrows; irow++) {
const Int_t off1 = irow*fNcols;
Int_t off2 = 0;
for (Int_t icol = 0; icol < fNcols; icol++)
data[off2+irow] = elem[off1+icol];
off2 += fNrows;
}
}
else
memcpy(data,elem,fNelems*sizeof(Element));
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::InsertRow(Int_t rown,Int_t coln,const Element *v,Int_t n)
{
const Int_t arown = rown-fRowLwb;
const Int_t acoln = coln-fColLwb;
const Int_t nr = (n > 0) ? n : fNcols;
if (gMatrixCheck) {
if (arown >= fNrows || arown < 0) {
Error("InsertRow","row %d out of matrix range",rown);
return *this;
}
if (acoln >= fNcols || acoln < 0) {
Error("InsertRow","column %d out of matrix range",coln);
return *this;
}
if (acoln+nr >= fNcols || nr < 0) {
Error("InsertRow","row length %d out of range",nr);
return *this;
}
}
const Int_t off = arown*fNcols+acoln;
Element * const elem = GetMatrixArray()+off;
memcpy(elem,v,nr*sizeof(Element));
return *this;
}
template<class Element>
void TMatrixTBase<Element>::ExtractRow(Int_t rown,Int_t coln,Element *v,Int_t n) const
{
const Int_t arown = rown-fRowLwb;
const Int_t acoln = coln-fColLwb;
const Int_t nr = (n > 0) ? n : fNcols;
if (gMatrixCheck) {
if (arown >= fNrows || arown < 0) {
Error("ExtractRow","row %d out of matrix range",rown);
return;
}
if (acoln >= fNcols || acoln < 0) {
Error("ExtractRow","column %d out of matrix range",coln);
return;
}
if (acoln+n >= fNcols || nr < 0) {
Error("ExtractRow","row length %d out of range",nr);
return;
}
}
const Int_t off = arown*fNcols+acoln;
const Element * const elem = GetMatrixArray()+off;
memcpy(v,elem,nr*sizeof(Element));
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::Shift(Int_t row_shift,Int_t col_shift)
{
fRowLwb += row_shift;
fColLwb += col_shift;
return *this;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::Zero()
{
R__ASSERT(IsValid());
memset(this->GetMatrixArray(),0,fNelems*sizeof(Element));
return *this;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::Abs()
{
R__ASSERT(IsValid());
Element *ep = this->GetMatrixArray();
const Element * const fp = ep+fNelems;
while (ep < fp) {
*ep = TMath::Abs(*ep);
ep++;
}
return *this;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::Sqr()
{
R__ASSERT(IsValid());
Element *ep = this->GetMatrixArray();
const Element * const fp = ep+fNelems;
while (ep < fp) {
*ep = (*ep) * (*ep);
ep++;
}
return *this;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::Sqrt()
{
R__ASSERT(IsValid());
Element *ep = this->GetMatrixArray();
const Element * const fp = ep+fNelems;
while (ep < fp) {
*ep = TMath::Sqrt(*ep);
ep++;
}
return *this;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::UnitMatrix()
{
R__ASSERT(IsValid());
Element *ep = this->GetMatrixArray();
memset(ep,0,fNelems*sizeof(Element));
for (Int_t i = fRowLwb; i <= fRowLwb+fNrows-1; i++)
for (Int_t j = fColLwb; j <= fColLwb+fNcols-1; j++)
*ep++ = (i==j ? 1.0 : 0.0);
return *this;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::NormByDiag(const TVectorT<Element> &v,Option_t *option)
{
R__ASSERT(IsValid());
R__ASSERT(v.IsValid());
if (gMatrixCheck) {
const Int_t nMax = TMath::Max(fNrows,fNcols);
if (v.GetNoElements() < nMax) {
Error("NormByDiag","vector shorter than matrix diagonal");
return *this;
}
}
TString opt(option);
opt.ToUpper();
const Int_t divide = (opt.Contains("D")) ? 1 : 0;
const Element *pV = v.GetMatrixArray();
Element *mp = this->GetMatrixArray();
if (divide) {
for (Int_t irow = 0; irow < fNrows; irow++) {
if (pV[irow] != 0.0) {
for (Int_t icol = 0; icol < fNcols; icol++) {
if (pV[icol] != 0.0) {
const Element val = TMath::Sqrt(TMath::Abs(pV[irow]*pV[icol]));
*mp++ /= val;
} else {
Error("NormbyDiag","vector element %d is zero",icol);
mp++;
}
}
} else {
Error("NormbyDiag","vector element %d is zero",irow);
mp += fNcols;
}
}
} else {
for (Int_t irow = 0; irow < fNrows; irow++) {
for (Int_t icol = 0; icol < fNcols; icol++) {
const Element val = TMath::Sqrt(TMath::Abs(pV[irow]*pV[icol]));
*mp++ *= val;
}
}
}
return *this;
}
template<class Element>
Element TMatrixTBase<Element>::RowNorm() const
{
R__ASSERT(IsValid());
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
Element norm = 0;
while (ep < fp) {
Element sum = 0;
for (Int_t j = 0; j < fNcols; j++)
sum += TMath::Abs(*ep++);
norm = TMath::Max(norm,sum);
}
R__ASSERT(ep == fp);
return norm;
}
template<class Element>
Element TMatrixTBase<Element>::ColNorm() const
{
R__ASSERT(IsValid());
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNcols;
Element norm = 0;
while (ep < fp) {
Element sum = 0;
for (Int_t i = 0; i < fNrows; i++,ep += fNcols)
sum += TMath::Abs(*ep);
ep -= fNelems-1;
norm = TMath::Max(norm,sum);
}
R__ASSERT(ep == fp);
return norm;
}
template<class Element>
Element TMatrixTBase<Element>::E2Norm() const
{
R__ASSERT(IsValid());
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
Element sum = 0;
for ( ; ep < fp; ep++)
sum += (*ep) * (*ep);
return sum;
}
template<class Element>
Int_t TMatrixTBase<Element>::NonZeros() const
{
R__ASSERT(IsValid());
Int_t nr_nonzeros = 0;
const Element *ep = this->GetMatrixArray();
const Element * const fp = ep+fNelems;
while (ep < fp)
if (*ep++ != 0.0) nr_nonzeros++;
return nr_nonzeros;
}
template<class Element>
Element TMatrixTBase<Element>::Sum() const
{
R__ASSERT(IsValid());
Element sum = 0.0;
const Element *ep = this->GetMatrixArray();
const Element * const fp = ep+fNelems;
while (ep < fp)
sum += *ep++;
return sum;
}
template<class Element>
Element TMatrixTBase<Element>::Min() const
{
R__ASSERT(IsValid());
const Element * const ep = this->GetMatrixArray();
const Int_t index = TMath::LocMin(fNelems,ep);
return ep[index];
}
template<class Element>
Element TMatrixTBase<Element>::Max() const
{
R__ASSERT(IsValid());
const Element * const ep = this->GetMatrixArray();
const Int_t index = TMath::LocMax(fNelems,ep);
return ep[index];
}
template<class Element>
void TMatrixTBase<Element>::Draw(Option_t *option)
{
gROOT->ProcessLine(Form("THistPainter::PaintSpecialObjects((TObject*)0x%x,\"%s\");",this,option));
}
template<class Element>
void TMatrixTBase<Element>::Print(Option_t *option) const
{
if (!IsValid()) {
Error("Print","Matrix is invalid");
return;
}
const char *format = "%11.4g ";
if (option) {
const char *f = strstr(option,"f=");
if (f) format = f+2;
}
Int_t i;
char topbar[100];
sprintf(topbar,format,123.456789);
Int_t nch = strlen(topbar)+1;
if (nch > 19) nch=19;
char ftopbar[20];
for (i=0;i<nch;i++) ftopbar[i] = ' ';
Int_t nk = 1 + Int_t(TMath::Log10(fNcols));
sprintf(ftopbar+nch/2,"%s%dd","%",nk);
Int_t nch2 = strlen(ftopbar);
for (i=nch2;i<nch;i++) ftopbar[i] = ' ';
ftopbar[nch] = '|';
ftopbar[nch+1] = 0;
printf("\n%dx%d matrix is as follows",fNrows,fNcols);
Int_t cols_per_sheet = 5;
if (nch <=8) cols_per_sheet =10;
const Int_t ncols = fNcols;
const Int_t nrows = fNrows;
const Int_t collwb = fColLwb;
const Int_t rowlwb = fRowLwb;
nk = 5+nch*TMath::Min(cols_per_sheet,fNcols);
for (i=0;i<nk;i++) topbar[i] = '-';
topbar[nk] = 0;
for (Int_t sheet_counter = 1; sheet_counter <= ncols; sheet_counter += cols_per_sheet) {
printf("\n\n |");
for (Int_t j = sheet_counter; j < sheet_counter+cols_per_sheet && j <= ncols; j++)
printf(ftopbar,j+collwb-1);
printf("\n%s\n",topbar);
for (Int_t i = 1; i <= nrows; i++) {
printf("%4d |",i+rowlwb-1);
for (Int_t j = sheet_counter; j < sheet_counter+cols_per_sheet && j <= ncols; j++)
printf(format,(*this)(i+rowlwb-1,j+collwb-1));
printf("\n");
}
}
printf("\n");
}
template<class Element>
Bool_t TMatrixTBase<Element>::operator==(Element val) const
{
R__ASSERT(IsValid());
if (val == 0. && fNelems == 0)
return kTRUE;
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
for (; ep < fp; ep++)
if (!(*ep == val))
return kFALSE;
return kTRUE;
}
template<class Element>
Bool_t TMatrixTBase<Element>::operator!=(Element val) const
{
R__ASSERT(IsValid());
if (val == 0. && fNelems == 0)
return kFALSE;
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
for (; ep < fp; ep++)
if (!(*ep != val))
return kFALSE;
return kTRUE;
}
template<class Element>
Bool_t TMatrixTBase<Element>::operator<(Element val) const
{
R__ASSERT(IsValid());
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
for (; ep < fp; ep++)
if (!(*ep < val))
return kFALSE;
return kTRUE;
}
template<class Element>
Bool_t TMatrixTBase<Element>::operator<=(Element val) const
{
R__ASSERT(IsValid());
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
for (; ep < fp; ep++)
if (!(*ep <= val))
return kFALSE;
return kTRUE;
}
template<class Element>
Bool_t TMatrixTBase<Element>::operator>(Element val) const
{
R__ASSERT(IsValid());
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
for (; ep < fp; ep++)
if (!(*ep > val))
return kFALSE;
return kTRUE;
}
template<class Element>
Bool_t TMatrixTBase<Element>::operator>=(Element val) const
{
R__ASSERT(IsValid());
const Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
for (; ep < fp; ep++)
if (!(*ep >= val))
return kFALSE;
return kTRUE;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::Apply(const TElementActionT<Element> &action)
{
R__ASSERT(IsValid());
Element *ep = this->GetMatrixArray();
const Element * const ep_last = ep+fNelems;
while (ep < ep_last)
action.Operation(*ep++);
return *this;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::Apply(const TElementPosActionT<Element> &action)
{
R__ASSERT(IsValid());
Element *ep = this->GetMatrixArray();
for (action.fI = fRowLwb; action.fI < fRowLwb+fNrows; action.fI++)
for (action.fJ = fColLwb; action.fJ < fColLwb+fNcols; action.fJ++)
action.Operation(*ep++);
R__ASSERT(ep == this->GetMatrixArray()+fNelems);
return *this;
}
template<class Element>
TMatrixTBase<Element> &TMatrixTBase<Element>::Randomize(Element alpha,Element beta,Double_t &seed)
{
R__ASSERT(IsValid());
const Element scale = beta-alpha;
const Element shift = alpha/scale;
Element * ep = GetMatrixArray();
const Element * const fp = ep+fNelems;
while (ep < fp)
*ep++ = scale*(Drand(seed)+shift);
return *this;
}
template<class Element>
Bool_t operator==(const TMatrixTBase<Element> &m1,const TMatrixTBase<Element> &m2)
{
if (!AreCompatible(m1,m2)) return kFALSE;
return (memcmp(m1.GetMatrixArray(),m2.GetMatrixArray(),
m1.GetNoElements()*sizeof(Element)) == 0);
}
template<class Element>
Element E2Norm(const TMatrixTBase<Element> &m1,const TMatrixTBase<Element> &m2)
{
if (gMatrixCheck && !AreCompatible(m1,m2)) {
::Error("E2Norm","matrices not compatible");
return -1.0;
}
const Element * mp1 = m1.GetMatrixArray();
const Element * mp2 = m2.GetMatrixArray();
const Element * const fmp1 = mp1+m1.GetNoElements();
Element sum = 0.0;
for (; mp1 < fmp1; mp1++, mp2++)
sum += (*mp1 - *mp2)*(*mp1 - *mp2);
return sum;
}
template<class Element1,class Element2>
Bool_t AreCompatible(const TMatrixTBase<Element1> &m1,const TMatrixTBase<Element2> &m2,Int_t verbose)
{
if (!m1.IsValid()) {
if (verbose)
::Error("AreCompatible", "matrix 1 not valid");
return kFALSE;
}
if (!m2.IsValid()) {
if (verbose)
::Error("AreCompatible", "matrix 2 not valid");
return kFALSE;
}
if (m1.GetNrows() != m2.GetNrows() || m1.GetNcols() != m2.GetNcols() ||
m1.GetRowLwb() != m2.GetRowLwb() || m1.GetColLwb() != m2.GetColLwb()) {
if (verbose)
::Error("AreCompatible", "matrices 1 and 2 not compatible");
return kFALSE;
}
return kTRUE;
}
template<class Element>
void Compare(const TMatrixTBase<Element> &m1,const TMatrixTBase<Element> &m2)
{
if (!AreCompatible(m1,m2)) {
Error("Compare(const TMatrixTBase<Element> &,const TMatrixTBase<Element> &)","matrices are incompatible");
return;
}
printf("\n\nComparison of two TMatrices:\n");
Element norm1 = 0;
Element norm2 = 0;
Element ndiff = 0;
Int_t imax = 0;
Int_t jmax = 0;
Element difmax = -1;
for (Int_t i = m1.GetRowLwb(); i <= m1.GetRowUpb(); i++) {
for (Int_t j = m1.GetColLwb(); j < m1.GetColUpb(); j++) {
const Element mv1 = m1(i,j);
const Element mv2 = m2(i,j);
const Element diff = TMath::Abs(mv1-mv2);
if (diff > difmax) {
difmax = diff;
imax = i;
jmax = j;
}
norm1 += TMath::Abs(mv1);
norm2 += TMath::Abs(mv2);
ndiff += TMath::Abs(diff);
}
}
printf("\nMaximal discrepancy \t\t%g", difmax);
printf("\n occured at the point\t\t(%d,%d)",imax,jmax);
const Element mv1 = m1(imax,jmax);
const Element mv2 = m2(imax,jmax);
printf("\n Matrix 1 element is \t\t%g", mv1);
printf("\n Matrix 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||Matrix 1|| \t\t\t%g", norm1);
printf("\n||Matrix 2|| \t\t\t%g", norm2);
printf("\n||Matrix1-Matrix2||\t\t\t\t%g", ndiff);
printf("\n||Matrix1-Matrix2||/sqrt(||Matrix1|| ||Matrix2||)\t%g\n\n",
ndiff/TMath::Max(TMath::Sqrt(norm1*norm2),1e-7));
}
template<class Element>
Bool_t VerifyMatrixValue(const TMatrixTBase<Element> &m,Element val,Int_t verbose,Element maxDevAllow)
{
R__ASSERT(m.IsValid());
if (m == 0)
return kTRUE;
Int_t imax = 0;
Int_t jmax = 0;
Element maxDevObs = 0;
if (TMath::Abs(maxDevAllow) <= 0.0)
maxDevAllow = std::numeric_limits<Element>::epsilon();
for (Int_t i = m.GetRowLwb(); i <= m.GetRowUpb(); i++) {
for (Int_t j = m.GetColLwb(); j <= m.GetColUpb(); j++) {
const Element dev = TMath::Abs(m(i,j)-val);
if (dev > maxDevObs) {
imax = i;
jmax = j;
maxDevObs = dev;
}
}
}
if (maxDevObs == 0)
return kTRUE;
if (verbose) {
printf("Largest dev for (%d,%d); dev = |%g - %g| = %g\n",imax,jmax,m(imax,jmax),val,maxDevObs);
if(maxDevObs > maxDevAllow)
Error("VerifyElementValue","Deviation > %g\n",maxDevAllow);
}
if(maxDevObs > maxDevAllow)
return kFALSE;
return kTRUE;
}
template<class Element>
Bool_t VerifyMatrixIdentity(const TMatrixTBase<Element> &m1,const TMatrixTBase<Element> &m2,Int_t verbose,
Element maxDevAllow)
{
if (!AreCompatible(m1,m2,verbose))
return kFALSE;
if (m1 == 0 && m2 == 0)
return kTRUE;
Int_t imax = 0;
Int_t jmax = 0;
Element maxDevObs = 0;
if (TMath::Abs(maxDevAllow) <= 0.0)
maxDevAllow = std::numeric_limits<Element>::epsilon();
for (Int_t i = m1.GetRowLwb(); i <= m1.GetRowUpb(); i++) {
for (Int_t j = m1.GetColLwb(); j <= m1.GetColUpb(); j++) {
const Element dev = TMath::Abs(m1(i,j)-m2(i,j));
if (dev > maxDevObs) {
imax = i;
jmax = j;
maxDevObs = dev;
}
}
}
if (maxDevObs == 0)
return kTRUE;
if (verbose) {
printf("Largest dev for (%d,%d); dev = |%g - %g| = %g\n",
imax,jmax,m1(imax,jmax),m2(imax,jmax),maxDevObs);
if (maxDevObs > maxDevAllow)
Error("VerifyMatrixValue","Deviation > %g\n",maxDevAllow);
}
if (maxDevObs > maxDevAllow)
return kFALSE;
return kTRUE;
}
template<class Element>
void TMatrixTBase<Element>::Streamer(TBuffer &R__b)
{
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) {
TMatrixTBase<Element>::Class()->ReadBuffer(R__b,this,R__v,R__s,R__c);
} else {
Error("TMatrixTBase<Element>::Streamer","Unknown version number: %d",R__v);
R__ASSERT(0);
}
if (R__v < 4) MakeValid();
} else {
TMatrixTBase<Element>::Class()->WriteBuffer(R__b,this);
}
}
template class TMatrixTBase<Float_t>;
template Bool_t operator== <Float_t>(const TMatrixFBase &m1,const TMatrixFBase &m2);
template Float_t E2Norm <Float_t>(const TMatrixFBase &m1,const TMatrixFBase &m2);
template Bool_t AreCompatible<Float_t,Float_t>
(const TMatrixFBase &m1,const TMatrixFBase &m2,Int_t verbose);
template Bool_t AreCompatible<Float_t,Double_t>
(const TMatrixFBase &m1,const TMatrixDBase &m2,Int_t verbose);
template void Compare <Float_t>(const TMatrixFBase &m1,const TMatrixFBase &m2);
template Bool_t VerifyMatrixValue <Float_t>(const TMatrixFBase &m,Float_t val,Int_t verbose,Float_t maxDevAllow);
template Bool_t VerifyMatrixValue <Float_t>(const TMatrixFBase &m,Float_t val);
template Bool_t VerifyMatrixIdentity<Float_t>(const TMatrixFBase &m1,const TMatrixFBase &m2,
Int_t verbose,Float_t maxDevAllowN);
template Bool_t VerifyMatrixIdentity<Float_t>(const TMatrixFBase &m1,const TMatrixFBase &m2);
template class TMatrixTBase<Double_t>;
template Bool_t operator== <Double_t>(const TMatrixDBase &m1,const TMatrixDBase &m2);
template Double_t E2Norm <Double_t>(const TMatrixDBase &m1,const TMatrixDBase &m2);
template Bool_t AreCompatible<Double_t,Double_t>
(const TMatrixDBase &m1,const TMatrixDBase &m2,Int_t verbose);
template Bool_t AreCompatible<Double_t,Float_t>
(const TMatrixDBase &m1,const TMatrixFBase &m2,Int_t verbose);
template void Compare <Double_t>(const TMatrixDBase &m1,const TMatrixDBase &m2);
template Bool_t VerifyMatrixValue <Double_t>(const TMatrixDBase &m,Double_t val,Int_t verbose,Double_t maxDevAllow);
template Bool_t VerifyMatrixValue <Double_t>(const TMatrixDBase &m,Double_t val);
template Bool_t VerifyMatrixIdentity<Double_t>(const TMatrixDBase &m1,const TMatrixDBase &m2,
Int_t verbose,Double_t maxDevAllow);
template Bool_t VerifyMatrixIdentity<Double_t>(const TMatrixDBase &m1,const TMatrixDBase &m2);
ROOT page - Class index - Class Hierarchy - Top of the page
This page has been automatically generated. If you have any comments or suggestions about the page layout send a mail to ROOT support, or contact the developers with any questions or problems regarding ROOT.