include/root/TDecompSVD.h

0001 // @(#)root/matrix:$Id$
0002 // Authors: Fons Rademakers, Eddy Offermann   Dec 2003
0003
0004 /*************************************************************************
0005  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
0006  * All rights reserved.                                                  *
0007  *                                                                       *
0008  * For the licensing terms see $ROOTSYS/LICENSE.                         *
0009  * For the list of contributors see $ROOTSYS/README/CREDITS.             *
0010  *************************************************************************/
0011
0012 #ifndef ROOT_TDecompSVD
0013 #define ROOT_TDecompSVD
0014
0015 ///////////////////////////////////////////////////////////////////////////
0016 //                                                                       //
0017 // Single Value Decomposition class                                      //
0018 //                                                                       //
0019 ///////////////////////////////////////////////////////////////////////////
0020
0021 #include "TDecompBase.h"
0022
0023 class TDecompSVD : public TDecompBase
0024 {
0025 protected :
0026
0027    //  A = fU fSig fV^T
0028    TMatrixD fU;    // orthogonal matrix
0029    TMatrixD fV;    // orthogonal matrix
0030    TVectorD fSig;  // diagonal of diagonal matrix
0031
0032    static Bool_t Bidiagonalize(TMatrixD &v,TMatrixD &u,TVectorD &sDiag,TVectorD &oDiag);
0033    static Bool_t Diagonalize  (TMatrixD &v,TMatrixD &u,TVectorD &sDiag,TVectorD &oDiag);
0034    static void   Diag_1       (TMatrixD &v,TVectorD &sDiag,TVectorD &oDiag,Int_t k);
0035    static void   Diag_2       (TVectorD &sDiag,TVectorD &oDiag,Int_t k,Int_t l);
0036    static void   Diag_3       (TMatrixD &v,TMatrixD &u,TVectorD &sDiag,TVectorD &oDiag,Int_t k,Int_t l);
0037    static void   SortSingular (TMatrixD &v,TMatrixD &u,TVectorD &sDiag);
0038
0039    const TMatrixDBase &GetDecompMatrix() const override { return fU; }
0040
0041 public :
0042
0043    enum {kWorkMax = 100}; // size of work array
0044
0045    TDecompSVD(): fU(), fV(), fSig() {}
0046    TDecompSVD(Int_t nrows,Int_t ncols);
0047    TDecompSVD(Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb);
0048    TDecompSVD(const TMatrixD &m,Double_t tol = 0.0);
0049    TDecompSVD(const TDecompSVD &another);
0050    ~TDecompSVD() override {}
0051
0052            const TMatrixD  GetMatrix ();
0053    Int_t     GetNrows  () const override;
0054    Int_t     GetNcols  () const override;
0055            const TMatrixD &GetU      ()       { if ( !TestBit(kDecomposed) ) Decompose();
0056                                                  return fU; }
0057            const TMatrixD &GetV      ()       { if ( !TestBit(kDecomposed) ) Decompose();
0058                                                 return fV; }
0059            const TVectorD &GetSig    ()       { if ( !TestBit(kDecomposed) ) Decompose();
0060                                                 return fSig; }
0061
0062    virtual       void      SetMatrix (const TMatrixD &a);
0063
0064    Bool_t   Decompose  () override;
0065    Bool_t   Solve      (      TVectorD &b) override;
0066    TVectorD Solve      (const TVectorD& b,Bool_t &ok) override { TVectorD x = b; ok = Solve(x);
0067                                                                 const Int_t rowLwb = GetRowLwb();
0068                                                                 x.ResizeTo(rowLwb,rowLwb+GetNcols()-1);
0069                                                                 return x; }
0070    Bool_t   Solve      (      TMatrixDColumn &b) override;
0071    Bool_t   TransSolve (      TVectorD &b) override;
0072    TVectorD TransSolve (const TVectorD& b,Bool_t &ok) override { TVectorD x = b; ok = TransSolve(x);
0073                                                                 const Int_t rowLwb = GetRowLwb();
0074                                                                 x.ResizeTo(rowLwb,rowLwb+GetNcols()-1);
0075                                                                 return x; }
0076    Bool_t   TransSolve (      TMatrixDColumn &b) override;
0077    Double_t Condition  () override;
0078    void     Det        (Double_t &d1,Double_t &d2) override;
0079
0080            Bool_t   Invert     (TMatrixD &inv);
0081            TMatrixD Invert     (Bool_t &status);
0082            TMatrixD Invert     () {Bool_t status; return Invert(status); }
0083
0084    void Print(Option_t *opt ="") const override; // *MENU*
0085
0086    TDecompSVD &operator= (const TDecompSVD &source);
0087
0088    ClassDefOverride(TDecompSVD,1) // Matrix Decompositition SVD
0089 };
0090
0091 #endif