From e2fdb0f4eac536de199baff2fa365f4a9b4bbe78 Mon Sep 17 00:00:00 2001 From: Pere Mato <pere.mato@cern.ch> Date: Thu, 27 Aug 2015 00:05:00 +0200 Subject: [PATCH] Doxygen continuation for hist --- .../doxygen/images/multidimfit_img86.gif | Bin 0 -> 3943 bytes hist/doc/hist.md | 4 + hist/hist/inc/THn.h | 57 +- hist/hist/inc/THnSparse.h | 58 +- hist/hist/inc/TLimit.h | 11 - hist/hist/src/TH2Poly.cxx | 8 +- hist/hist/src/TH3.cxx | 17 +- hist/hist/src/THLimitsFinder.cxx | 1 + hist/hist/src/THStack.cxx | 39 +- hist/hist/src/THn.cxx | 1 + hist/hist/src/THnBase.cxx | 31 +- hist/hist/src/THnSparse.cxx | 2 + hist/hist/src/TKDE.cxx | 18 +- hist/hist/src/TLimit.cxx | 111 +- hist/hist/src/TLimitDataSource.cxx | 8 +- hist/hist/src/TMultiDimFit.cxx | 1948 +++-------------- hist/hist/src/TMultiGraph.cxx | 35 +- hist/hist/src/TPolyMarker.cxx | 14 +- hist/hist/src/TPrincipal.cxx | 620 +----- hist/hist/src/TProfile.cxx | 564 +++-- hist/hist/src/TProfile2D.cxx | 253 +-- hist/hist/src/TProfile3D.cxx | 202 +- hist/hist/src/TSVDUnfold.cxx | 9 +- hist/hist/src/TSpline.cxx | 12 +- hist/hist/src/TUnfold.cxx | 453 ++-- hist/hist/src/TUnfoldBinning.cxx | 59 +- hist/hist/src/TUnfoldDensity.cxx | 202 +- hist/hist/src/TUnfoldSys.cxx | 193 +- hist/hist/src/TVirtualFitter.cxx | 7 +- hist/hist/src/TVirtualGraphPainter.cxx | 9 +- hist/hist/src/TVirtualHistPainter.cxx | 8 +- 31 files changed, 1500 insertions(+), 3454 deletions(-) create mode 100644 documentation/doxygen/images/multidimfit_img86.gif create mode 100644 hist/doc/hist.md diff --git a/documentation/doxygen/images/multidimfit_img86.gif b/documentation/doxygen/images/multidimfit_img86.gif new file mode 100644 index 0000000000000000000000000000000000000000..b0353b4cda2aff90770ae8f0d5816de332d94d4b GIT binary patch literal 3943 zcma)%=RX??!-f+QNstgTMvcUtv1-qlF`7E4Q6;sjRaC7Qu|;jM+gh<VrRXqQQuLHk zDMxFxMYTptk5zqs=lu)b`@{9=`f%MAR_0pT#93eqpce4o+S}V}YHBJfDoRUBQ&Lh$ zBvM2~go}%dqobp>wY8z4;m*#EwzjsKnwp}bqO7c}sHo`F)RcgL06#xJFE1|=iG;)9 z5C{Ye1~)Y|(P%U(m70-}k)EFJ>+9?3>1l6oPaqKf{P|;PX=!R|YG7cnyu6G?qlJZq zQ79BACno}d7#tjAXJ>!%<OviCWnyCbzpnom#E(Hh0009}_|N|DlM&1LP+qX5f?&=B zEtGy!r^1HLhHrIkXS9S*evckp>4^NM2S7;Pq3nsP0zlu{2%)3&qo;9~G-~E=n%xd@ zkqk3K0AMkhHV|n`X-F2eu&CH62P4g8o?lW~RqYBXz$wd!3L@%EYHO@~{jx*=y~!Q@ zw(U6odca?8tpjbnqokh4W0O;t03>D=rs=tP&%dMDmH2btefVg?APVv=uWW4Et$BX> zvb7rx6LHvfVdop(qy9LtJ$r5Yk7R#wyJ$uH9AWB-JqKPliOOT*fS$pR+5&W<-Z1u7 zO$~S(#iv3;#GdD5nOE~v&!V~AOrULgeOThOl4P2H#|t&-b1qqSUmch%b(60GqsKb! z73N~CeCR!k`B7zO?`fT*XV!n%P+8!zzU=+(fXxtqw-B@K3J-9h-hMKHbIZ0|>Tp^c z$tR;XQ7HnaYxAQE?QcdR33@oa5_Yc}%EI3$|7dZ1ZYe(WZWH!Z>mH9Zx^s>e^x4AJ zTimvSPpSzlH1!XiYYEwWXE+D#@@Af12uvCzX|gHd9|!xsJd`wS3)Gf*1LNry3=zv? zZl0?6Y49}GGz{>zIcfaiZ?1hpcSP9|C2&l5^htJM{F<J7G+=x8EOFQJ{7mYwbYVHZ z0`6ay0Qwp#A_J7Nv`Shr+7uJzvP)lNLp>y~vKgNf8CLw#b8CP!K2>S6RH5*rk+f?l zlMm30<11@6(v_0Q5J9s4I#lbmrDvA*4-_p2KL@QyOZ*3n&niLXGUZvJ-dkUgmkNLv z*cw=QGU2ktQaB}~deSxVp+plfL1(xYsSKqJm&{Z9!9d<ed0$LJoscHb2;X99UeV;C z>lM=yZ$5Cv%NYYeNhkg?*K#iAHE^*Qz44>XwXEjzi87+<>xc=gV0(69Y!UJSm)s3C zNmS*FibAfhpPSTr@cHg~u>IjG6NN<?YOfP+yVaH^w+$VwrhcolLCJ4-^@B)ouIJCg zE{B;|l6n1LsJ-T{Han7O&p28rsP6n}nPhfneA}pzah9=e-ct+YNKpSlMi-FnEVil< zr0|Ha!+HAa6gk9HP!%xDOcLnO5Uk{MmDbfzFdZWj^*UX-G+R_OxSTz)H_&dvKPEUI z*f`VS%{J5HQA1bTrsQW7Qm2#GTEr&#c2|E8_5Qs-{`d5zg1YVnWuJqYYYO_TGdA~t zJpslRwyLi>7+2ZuE!G7da@lRy+%?CA{np?mrm6ODE5`taAhiH;kGv#Je0!05b=q|4 z_PD<0-Mf=oH6L80V=vAM2*pLd%}~CU)yvn=)4Faj;M@!5yN}n#fP_anvGp;zWT`7t z#@Ci=*d86LK5zpYjXx#J>YQ&mIol<R7qqpdh!><gCqwckSC6xXw|C&XdWn9ZJvYXV z4?0mu`S-z`sQB<b3+mpj!$@~SNv_2^A9En1&K*!@nMZz28n92<DFSlpHU>z@o^OP5 z*~aX-<eqd%serMu$2W6RxYZHxT4cI&qe5@#$?*%|LqQ<0lsKd&hE0g2d<a3T^aurD z`FBW6%yME(%^|_CpBQe8QDveMB?l?)od`WT9JxKlc~c%>3|9^r#bAekO8<y!MEfU; zvRqEb@{9t}>8634o%}X8*LHLwLFkdzOUGX=_(p-L@^Kj%C#RLj>kUKJl*^(u<@oLr z!b{!S86hL@bp7)+hP00Fi|m>V5Nqh1E}aLNLtye&69$-@x)Ui~Xo+h)vg7FdR`?J% z6b*Y}8F=$DEb2#f6Ki9RHCkL6YdkN&{mIJa>5ROcKc>n~mdEM3u~n1$x+cg<K@!z@ z_0gH7QmIC+`^iLce}jYr6-!r1<$^sYOn&r5+xUj6=5(g7Jt(6$`8lnk%cu^>NAui< zPIc5xK>bCt_t+4pk;1hBRN_l|&BdpcN`=N7_)ipW^q9TVfdmV-{wUKebieEpXx!jb zg`j}C4mx+3>|_<pfbak~0$MKfaPtxZwbT`)Bgq=Ax#j>d&&=3Yo8~GZiFZTAvnAcX zWZIvX@fPI+Yo${1ty>dV^<;kB)o+&d6@;*uFeTN;&NgZ~b*4EwgsGcO<LwU#DmhWz zg3>EPf6<2$AO|f?#Nj59Y8^+HfJ3~>;fT$o#=AnQBfFoZgDXYcVG^Lxa684)&4YVE zeY(fr{lcg{q<5k8z_C4Q+QpMnl>lkj@lAU6_DV{Zv+CJ=v(xXFuC+4BbiuJM`>~sW zU*%prsUGX5>Vy5di=1B>{KC}3Z`tF^@4ajahqjqn-Sp61^6}@RJdU*pis*A=*+0-R z4F&{iLKZ2Je0@019mG9r#!xj`Pm;`SO?#63H8a$tF=Zm4f~5OZiaaoOAK{Uvd!E5~ zz5*yAB})i<CNO_uPQ;=9lKtYwdUY(D(36rLHM8`AyXEn--9u{3mfOc`SBs5~=gn?= z^mzuWmykHxKzC$ySK}P|DEuO3;kMQvB5pcoh`x5cWNx#@Rem%kW7nB+WonJ@$<Yh> z*QiWCn{~mZqnBFQ`+<ta<PFy$pxkn&fw#(a9$`<CiTP`Bb?mIf;Aqq}CNenp+b3t| zmAfEI5UbAc2WELWBp=kHP$1!&LSX6ScQ7+GFXgQYZWVOJ7ogpd$P^KCB~p2yyPOX2 z>z|Hdy{nY?5U!&qj;f9%1ts<p!~}3F`}&bmOf+u&R5MLMJ*E5VRrH^!=k&{Euqken zvLty=d2g3jj09zC19&PS$8(OBeVZ3bNR&_O2hCQgHgBz8>;n=mKKr%QCd>29gX5Gj z@p2kgkUQdiWj*<ZmPz-6*jx<Xs_n(8dCus>W-t8frGGy3_AbBe(3pY47~yMugU_;q z)4Y~<d>_+bk4o)&vMAg%b44X2ncb%5Dm>Y92?3~xQS<tE_DeHiN2p+eJ!y5U0U0Nl z!0KlC)9zo}FHjvnMEZuwTG|WTi;qQ>zJDih25d&2o6q>krp2(c5)$ftG;U=zL*l3X zwe<v$KTG{A@Z}w3IrE4?O-2Wjy@cE|FQ+Lpzz<JV8Q_*hlJ(DvCP=+JKj)6$>;kax zZI^{Zne4BlFdMOnkt%nj0T9u}p8!(ODqmdK@e%ov0Fdv~>M<^Zuy|G`lQ@!=(f~{- z1^vBnI&Kg`jvW3smWjtmw`Vz<m9F{l!S#l-&ng_L*Ms6py!OsE-lJ1V!f^*r4-&V8 zc0I`CDrVyEOWP_ZKhU59u;1X@Z_kXimg0f2{ErB`rz&UhO(Njo3c|kG%#z42En4To zco5|%X4%4nhix&OzsrjDt$y_M){L&>S0&+HFZAg11*Q|b-N6~T@!u4W5vE%Khq_gU zJOS&Zpf<ybGZN~2e*;%bU=qm%h}TC&ti@MyUaz2@gR3nh*r6QAk~nkhg$gXhVvlT~ zn!r9wbPNMG^E-(6i(*jZ{8V_js=o>m_-(^bG7<*+#pOgM6^#q0IFeY9NhRvnaYlh4 ztPX5QwBt(BcvaG+zcdycz(Wz5mXTQ7&FEe~c7x;?8_Q%WCsB$4ny+^{OG7yO$2dB0 zxItr`7|<^Z<V!=b0imfMhf;a^k}#F2gtb%-CUWW=XgfmT-!uftGCI;ZK6y5+xRM;v z2ZENUg-0TiH`7uXDStZCi-r<Zo!!?oG5JtL@kIKJx|*^-jJKL%HkP=|7guLbvEHZ5 zlQT{>aMcZI51e7<4dxcu<Y0e&@ADcyx-`$#LdU#Au(t8NajtteekQl8ru3)FdwH|= z5K=cfVxAjkwtHtXuV-4Xr4txQ;hVQ!`g-5rPa*hZ-G3ZHIl`<s!@{&EII*-%%L{e6 z07e1t#*x<NzU&aoxSz-z9@XqiVj&XmEKQ6oT{TgODDjg9gycmmvkNSD6C|DlhccAL zJ9BZ4xf$=mIHDA|-a|__L9&PJ4=pY~>&(L)W&o>Db4cSPX&$syh=c$*QQSufl;2sE zkvo%bO+aCfQwb&oW<Ee$bI}Xk0PS!u7Gs3rF&xebL`$pO!BX|K1g)YFU*<r9YK4D@ zg(IQV#yJ(|(hKf%=(&3**!=Sa`3-ZQ08XN85hJ(Abfk!QAedhT{dFXI3_{-R2j}gx zN4zh(V+He|1jW%rKhp|}OiOO^6vxX7QyL-cg<#4XwJZuE?-;MBg(;xGqF2GCFv;>p zjt6vbg$ue$mUT{orPj(S)e2rOz|veQ)!N9mj1IplrpSyd?+I1n{H)05RlYd_x%XaZ zP(ZQqJ!n`!=!Gmx5FJd{#!QwXrYYdRwZ&(mU?f@oSQn9b7{~iY@RAF98OAyu09vIB ztiv9T(m|gElsD)lh5g{-Qn6i_#rLg}p8?|g0!Y8c(qk)$8Vcff6!_1U81M$HeG3f6 zOF&r9Ibcm-Hdi$KM0jgUvdcxgFRI)?tx9dA{F$~d@&@dJoVw@$_$^%@^A#?M*OwmP zSP0Pj7==)jGf<9(i4T--Td}L@lxzf7ysfEea;Z=rtI5Px&SPlrjA(DYXjyxehORXO zm9({WnyPD+NnMqwacz@N?PPFG51lr$<?9*^mX*^dy2_kVOueExf8+tuo!GScbX0Cs z2MpA4U%}mWg^k){MzL%f2er|`_0fU|d3=p_wCSd!=3Q>TssbImA;7H@YXdn#a&%zY zVpvyALnNc2=(wSnz_%F6-H}R26cB|c(WZwQ0V|E8SSVp7B*LUACb%pJXGqP*>&0?5 z^Z{#wKWTI>vUR`0q)RpRjn%SZYESzx`m>@(CM}LlEkI(6;zW~aV(mq_mbB6OD9_gN zPt`thrhZ>}mIiRkf*idA+IL;ypXKV-qg&^qFEX02+b0^|Pe4XZx<3cl4+pgWeZ}+9 r71uVv{(C_CPZ^{Sug~lT2g~b2CZSy9D&0pN#MKTjmd+JU5J2~TRffl# literal 0 HcmV?d00001 diff --git a/hist/doc/hist.md b/hist/doc/hist.md new file mode 100644 index 00000000000..febaec2af4c --- /dev/null +++ b/hist/doc/hist.md @@ -0,0 +1,4 @@ +/** +\defgroup Hist Histogram Library +The Histogram library is documented in the class TH1. +*/ \ No newline at end of file diff --git a/hist/hist/inc/THn.h b/hist/hist/inc/THn.h index 6f76faf28df..5d9369403d5 100644 --- a/hist/hist/inc/THn.h +++ b/hist/hist/inc/THn.h @@ -206,34 +206,35 @@ protected: //______________________________________________________________________________ -// -// Templated implementation of the abstract base THn. -// All functionality and the interfaces to be used are in THn! -// -// THn does not know how to store any bin content itself. Instead, this -// is delegated to the derived, templated class: the template parameter decides -// what the format for the bin content is. The actual storage is delegated to -// TNDArrayT<T>. -// -// Typedefs exist for template parematers with ROOT's generic types: -// -// Templated name Typedef Bin content type -// THnT<Char_t> THnC Char_t -// THnT<Short_t> THnS Short_t -// THnT<Int_t> THnI Int_t -// THnT<Long_t> THnL Long_t -// THnT<Float_t> THnF Float_t -// THnT<Double_t> THnD Double_t -// -// We recommend to use THnC wherever possible, and to map its value space -// of 256 possible values to e.g. float values outside the class. This saves an -// enourmous amount of memory. Only if more than 256 values need to be -// distinguished should e.g. THnS or even THnF be chosen. -// -// Implementation detail: the derived, templated class is kept extremely small -// on purpose. That way the (templated thus inlined) uses of this class will -// only create a small amount of machine code, in contrast to e.g. STL. -//______________________________________________________________________________ +/** \class THnT + Templated implementation of the abstract base THn. + All functionality and the interfaces to be used are in THn! + + THn does not know how to store any bin content itself. Instead, this + is delegated to the derived, templated class: the template parameter decides + what the format for the bin content is. The actual storage is delegated to + TNDArrayT<T>. + + Typedefs exist for template parematers with ROOT's generic types: + + Templated name | Typedef | Bin content type + -----------------|---------------|-------------------- + THnT<Char_t> | THnC | Char_t + THnT<Short_t> | THnS | Short_t + THnT<Int_t> | THnI | Int_t + THnT<Long_t> | THnL | Long_t + THnT<Float_t> | THnF | Float_t + THnT<Double_t> | THnD | Double_t + + We recommend to use THnC wherever possible, and to map its value space + of 256 possible values to e.g. float values outside the class. This saves an + enourmous amount of memory. Only if more than 256 values need to be + distinguished should e.g. THnS or even THnF be chosen. + + Implementation detail: the derived, templated class is kept extremely small + on purpose. That way the (templated thus inlined) uses of this class will + only create a small amount of machine code, in contrast to e.g. STL. +*/ template <typename T> class THnT: public THn { diff --git a/hist/hist/inc/THnSparse.h b/hist/hist/inc/THnSparse.h index e2d6423b48b..c4a2b99eb86 100644 --- a/hist/hist/inc/THnSparse.h +++ b/hist/hist/inc/THnSparse.h @@ -184,34 +184,36 @@ class THnSparse: public THnBase { //______________________________________________________________________________ -// -// Templated implementation of the abstract base THnSparse. -// All functionality and the interfaces to be used are in THnSparse! -// -// THnSparse does not know how to store any bin content itself. Instead, this -// is delegated to the derived, templated class: the template parameter decides -// what the format for the bin content is. In fact it even defines the array -// itself; possible implementations probably derive from TArray. -// -// Typedefs exist for template parematers with ROOT's generic types: -// -// Templated name Typedef Bin content type -// THnSparseT<TArrayC> THnSparseC Char_t -// THnSparseT<TArrayS> THnSparseS Short_t -// THnSparseT<TArrayI> THnSparseI Int_t -// THnSparseT<TArrayL> THnSparseL Long_t -// THnSparseT<TArrayF> THnSparseF Float_t -// THnSparseT<TArrayD> THnSparseD Double_t -// -// We recommend to use THnSparseC wherever possible, and to map its value space -// of 256 possible values to e.g. float values outside the class. This saves an -// enourmous amount of memory. Only if more than 256 values need to be -// distinguished should e.g. THnSparseS or even THnSparseF be chosen. -// -// Implementation detail: the derived, templated class is kept extremely small -// on purpose. That way the (templated thus inlined) uses of this class will -// only create a small amount of machine code, in contrast to e.g. STL. -//______________________________________________________________________________ +/** \class THnSparseT + Templated implementation of the abstract base THnSparse. + All functionality and the interfaces to be used are in THnSparse! + + THnSparse does not know how to store any bin content itself. Instead, this + is delegated to the derived, templated class: the template parameter decides + what the format for the bin content is. In fact it even defines the array + itself; possible implementations probably derive from TArray. + + Typedefs exist for template parematers with ROOT's generic types: + + Templated name | Typedef | Bin content type + -----------------|--------------|-------------------- + THnSparseT<TArrayC> | THnSparseC | Char_t + THnSparseT<TArrayS> | THnSparseS | Short_t + THnSparseT<TArrayI> | THnSparseI | Int_t + THnSparseT<TArrayL> | THnSparseL | Long_t + THnSparseT<TArrayF> | THnSparseF | Float_t + THnSparseT<TArrayD> | THnSparseD | Double_t + + We recommend to use THnSparseC wherever possible, and to map its value space + of 256 possible values to e.g. float values outside the class. This saves an + enourmous amount of memory. Only if more than 256 values need to be + distinguished should e.g. THnSparseS or even THnSparseF be chosen. + + Implementation detail: the derived, templated class is kept extremely small + on purpose. That way the (templated thus inlined) uses of this class will + only create a small amount of machine code, in contrast to e.g. STL. +*/ + template <class CONT> class THnSparseT: public THnSparse { diff --git a/hist/hist/inc/TLimit.h b/hist/hist/inc/TLimit.h index 09ca9f1c630..a1cedaf07d3 100644 --- a/hist/hist/inc/TLimit.h +++ b/hist/hist/inc/TLimit.h @@ -17,17 +17,6 @@ class TArrayD; class TOrdCollection; class TH1; -//____________________________________________________________________ -// -// TLimit -// -// This class computes 95% Confidence Levels. -// -// Implemented by C. Delaere from the mclimit code written by Tom Junk. -// reference: HEP-EX/9902006 -// See http://cern.ch/thomasj/searchlimits/ecl.html for more details. -//____________________________________________________________________ - class TLimit { protected: static bool Fluctuate(TLimitDataSource * input, TLimitDataSource * output, bool init,TRandom *, bool stat=false); diff --git a/hist/hist/src/TH2Poly.cxx b/hist/hist/src/TH2Poly.cxx index 587798ee5b4..37764b33ae0 100644 --- a/hist/hist/src/TH2Poly.cxx +++ b/hist/hist/src/TH2Poly.cxx @@ -28,7 +28,8 @@ ClassImp(TH2Poly) -/** \class TH2Poly +/** \class TH2Poly + \ingroup Hist 2D Histogram with Polygonal Bins <h3>Overview</h3> @@ -1242,6 +1243,11 @@ void TH2Poly::SetFloat(Bool_t flag) } +/** \class TH2PolyBin +Helper class to represent a bin in the TH2Poly histogram +*/ + + //////////////////////////////////////////////////////////////////////////////// /// Default constructor. diff --git a/hist/hist/src/TH3.cxx b/hist/hist/src/TH3.cxx index 943b58d36d6..1ba4b398b88 100644 --- a/hist/hist/src/TH3.cxx +++ b/hist/hist/src/TH3.cxx @@ -26,13 +26,18 @@ ClassImp(TH3) -/// \class TH3C \brief tomato 3-D histogram with a bype per channel (see TH1 documentation) -/// \class TH3S \brief tomato 3-D histogram with a short per channel (see TH1 documentation) -/// \class TH3I \brief tomato 3-D histogram with a int per channel (see TH1 documentation)} -/// \class TH3F \brief tomato 3-D histogram with a float per channel (see TH1 documentation)} -/// \class TH3D \brief tomato 3-D histogram with a double per channel (see TH1 documentation)} +/** \addtogroup Hist +@{ +\class TH3C \brief tomato 3-D histogram with a bype per channel (see TH1 documentation) +\class TH3S \brief tomato 3-D histogram with a short per channel (see TH1 documentation) +\class TH3I \brief tomato 3-D histogram with a int per channel (see TH1 documentation)} +\class TH3F \brief tomato 3-D histogram with a float per channel (see TH1 documentation)} +\class TH3D \brief tomato 3-D histogram with a double per channel (see TH1 documentation)} +@} +*/ -/** \class TH3 +/** \class TH3 + \ingroup Hist The 3-D histogram classes derived from the 1-D histogram classes. All operations are supported (fill, fit). Drawing is currently restricted to one single option. diff --git a/hist/hist/src/THLimitsFinder.cxx b/hist/hist/src/THLimitsFinder.cxx index 3c463f3285e..b25bb3be520 100644 --- a/hist/hist/src/THLimitsFinder.cxx +++ b/hist/hist/src/THLimitsFinder.cxx @@ -30,6 +30,7 @@ THLimitsFinder *THLimitsFinder::fgLimitsFinder = 0; ClassImp(THLimitsFinder) /** \class THLimitsFinder + \ingroup Hist Class to find nice axis limits */ diff --git a/hist/hist/src/THStack.cxx b/hist/hist/src/THStack.cxx index 3e2c4872d7f..5852e37c1f7 100644 --- a/hist/hist/src/THStack.cxx +++ b/hist/hist/src/THStack.cxx @@ -28,31 +28,31 @@ ClassImp(THStack) //////////////////////////////////////////////////////////////////////////////// -/* Begin_Html -<center><h2>The Histogram stack class</h2></center> +/** \class THStack + \ingroup Hist +The Histogram stack class -A <t>THStack</tt> is a collection of <tt>TH1</tt> or <tt>TH2</tt> histograms. -Using <tt>THStack::Draw()</tt> the histogram collection is drawn in one go according +A THStack is a collection of TH1 or TH2 histograms. +Using THStack::Draw() the histogram collection is drawn in one go according to the drawing option. -<p> -<tt>THStack::Add()</tt> allows to add a new histogram to the list. -The <tt>THStack</tt> does not own the objects in the list. -<p> + +THStack::Add() allows to add a new histogram to the list. +The THStack does not own the objects in the list. + By default (if no option drawing option is specified), histograms will be paint -stacked on top of each other. <tt>TH2</tt are stacked as lego plots. -<p> -If option <tt>"nostack"</tt> is specified the histograms are not drawn on top -of each other but as they would if drawn using the option <tt>"same"</tt>. -<p> -If option <tt>"nostackb"</tt> is specified the histograms are drawn next to +stacked on top of each other. TH2 are stacked as lego plots. + +If option "nostack" is specified the histograms are not drawn on top +of each other but as they would if drawn using the option "same". + +If option "nostackb" is specified the histograms are drawn next to each other as bar charts. -<p> + In all cases The axis range is computed automatically along the X and Y axis in order to show the complete histogram collection. -<p> + Example; -End_Html Begin_Macro(source) { THStack *hs = new THStack("hs",""); @@ -78,18 +78,15 @@ Begin_Macro(source) return cs; } End_Macro -Begin_Html A more complex example: -End_Html Begin_Macro(source) ../../../tutorials/hist/hstack.C End_Macro -Begin_Html Note that picking is supported for all drawing modes. -End_Html */ +*/ //////////////////////////////////////////////////////////////////////////////// diff --git a/hist/hist/src/THn.cxx b/hist/hist/src/THn.cxx index f7b1367c784..ed352f202d8 100644 --- a/hist/hist/src/THn.cxx +++ b/hist/hist/src/THn.cxx @@ -122,6 +122,7 @@ namespace { /** \class THn + \ingroup Hist Multidimensional histogram. Use a THn if you really, really have to store more than three dimensions, diff --git a/hist/hist/src/THnBase.cxx b/hist/hist/src/THnBase.cxx index 5635dbe6d9b..b6398accfb3 100644 --- a/hist/hist/src/THnBase.cxx +++ b/hist/hist/src/THnBase.cxx @@ -32,6 +32,7 @@ /** \class THnBase + \ingroup Hist Multidimensional histogram base. Defines common functionality and interfaces for THn, THnSparse. */ @@ -1391,12 +1392,11 @@ void THnBase::Browse(TBrowser *b) -//______________________________________________________________________________ -// -// -// Iterator over THnBase bins; internal implementation. -// -//////////////////////////////////////////////////////////////////////////////// + +/** \class ROOT::THnBaseBinIter + Iterator over THnBase bins (internal implementation). +*/ + /// Destruct a bin iterator. ROOT::THnBaseBinIter::~THnBaseBinIter() { @@ -1404,12 +1404,10 @@ ROOT::THnBaseBinIter::~THnBaseBinIter() { } -//______________________________________________________________________________ -// -// -// Iterator over THnBase bins -// -//////////////////////////////////////////////////////////////////////////////// + +/** \class THnIter + Iterator over THnBase bins +*/ ClassImp(THnIter); @@ -1419,14 +1417,11 @@ THnIter::~THnIter() { } +/** \class ROOT::THnBaseBrowsable + TBrowser helper for THnBase. +*/ -//______________________________________________________________________________ -// -// TBrowser helper for THnBase. -// -//////////////////////////////////////////////////////////////////////////////// - ClassImp(ROOT::THnBaseBrowsable); //////////////////////////////////////////////////////////////////////////////// diff --git a/hist/hist/src/THnSparse.cxx b/hist/hist/src/THnSparse.cxx index dfa3d86dd0e..b5600a747f4 100644 --- a/hist/hist/src/THnSparse.cxx +++ b/hist/hist/src/THnSparse.cxx @@ -490,6 +490,8 @@ void THnSparseArrayChunk::Sumw2() /** \class THnSparse + \ingroup Hist + Efficient multidimensional histogram. Use a THnSparse instead of TH1 / TH2 / TH3 / array for histogramming when diff --git a/hist/hist/src/TKDE.cxx b/hist/hist/src/TKDE.cxx index b372bcd996b..559e379044c 100644 --- a/hist/hist/src/TKDE.cxx +++ b/hist/hist/src/TKDE.cxx @@ -9,15 +9,19 @@ * * **********************************************************************/ ////////////////////////////////////////////////////////////////////////////// -/* +/** \class TKDE + \ingroup Hist Kernel Density Estimation class. - The three main references are (1) "Scott DW, Multivariate Density Estimation. Theory, Practice and Visualization. New York: Wiley", - (2) "Jann Ben - ETH Zurich, Switzerland -, Univariate kernel density estimation document for KDENS: - Stata module for univariate kernel density estimation." - (3) "Hardle W, Muller M, Sperlich S, Werwatz A, Nonparametric and Semiparametric Models. Springer." - The algorithm is briefly described in (4) "Cranmer KS, Kernel Estimation in High-Energy + The three main references are: + 1. "Scott DW, Multivariate Density Estimation. Theory, Practice and Visualization. New York: Wiley", + 2. "Jann Ben - ETH Zurich, Switzerland -, Univariate kernel density estimation document for KDENS: + Stata module for univariate kernel density estimation." + 3. "Hardle W, Muller M, Sperlich S, Werwatz A, Nonparametric and Semiparametric Models. Springer." + 4. "Cranmer KS, Kernel Estimation in High-Energy Physics. Computer Physics Communications 136:198-207,2001" - e-Print Archive: hep ex/0011057. - A binned version is also implemented to address the performance issue due to its data size dependance. + + The algorithm is briefly described in (4). A binned version is also implemented to address the + performance issue due to its data size dependance. */ diff --git a/hist/hist/src/TLimit.cxx b/hist/hist/src/TLimit.cxx index 9374d6e02b1..45ae6b539e3 100644 --- a/hist/hist/src/TLimit.cxx +++ b/hist/hist/src/TLimit.cxx @@ -28,6 +28,58 @@ // /////////////////////////////////////////////////////////////////////////// +/** \class TLimit + \ingroup Hist + Algorithm to compute 95% C.L. limits using the Likelihood ratio + semi-bayesian method. + + Implemented by C. Delaere from the mclimit code written by Tom Junk [HEP-EX/9902006]. + See [http://cern.ch/thomasj/searchlimits/ecl.html](http://cern.ch/thomasj/searchlimits/ecl.html) for more details. + + It takes signal, background and data histograms wrapped in a + TLimitDataSource as input and runs a set of Monte Carlo experiments in + order to compute the limits. If needed, inputs are fluctuated according + to systematics. The output is a TConfidenceLevel. + + The class TLimitDataSource takes the signal, background and data histograms as well as different + systematics sources to form the TLimit input. + + The class TConfidenceLevel represents the final result of the TLimit algorithm. It is created just after the + time-consuming part and can be stored in a TFile for further processing. + It contains light methods to return CLs, CLb and other interesting + quantities. + + The actual algorithm... + + From an input (TLimitDataSource) it produces an output TConfidenceLevel. + For this, nmc Monte Carlo experiments are performed. + As usual, the larger this number, the longer the compute time, + but the better the result. + + Supposing that there is a plotfile.root file containing 3 histograms + (signal, background and data), you can imagine doing things like: + +~~~{.cpp} + TFile* infile=new TFile("plotfile.root","READ"); + infile->cd(); + TH1* sh=(TH1*)infile->Get("signal"); + TH1* bh=(TH1*)infile->Get("background"); + TH1* dh=(TH1*)infile->Get("data"); + TLimitDataSource* mydatasource = new TLimitDataSource(sh,bh,dh); + TConfidenceLevel *myconfidence = TLimit::ComputeLimit(mydatasource,50000); + std::cout << " CLs : " << myconfidence->CLs() << std::endl; + std::cout << " CLsb : " << myconfidence->CLsb() << std::endl; + std::cout << " CLb : " << myconfidence->CLb() << std::endl; + std::cout << "< CLs > : " << myconfidence->GetExpectedCLs_b() << std::endl; + std::cout << "< CLsb > : " << myconfidence->GetExpectedCLsb_b() << std::endl; + std::cout << "< CLb > : " << myconfidence->GetExpectedCLb_b() << std::endl; + delete myconfidence; + delete mydatasource; + infile->Close(); +~~~ + More information can still be found on [this page](http://cern.ch/aleph-proj-alphapp/doc/tlimit.html) + */ + #include "TLimit.h" #include "TArrayD.h" #include "TVectorD.h" @@ -52,65 +104,6 @@ TConfidenceLevel *TLimit::ComputeLimit(TLimitDataSource * data, Int_t nmc, bool stat, TRandom * generator) { - // class TLimit - // ------------ - // Algorithm to compute 95% C.L. limits using the Likelihood ratio - // semi-bayesian method. - // It takes signal, background and data histograms wrapped in a - // TLimitDataSource as input and runs a set of Monte Carlo experiments in - // order to compute the limits. If needed, inputs are fluctuated according - // to systematics. The output is a TConfidenceLevel. - // - // class TLimitDataSource - // ---------------------- - // - // Takes the signal, background and data histograms as well as different - // systematics sources to form the TLimit input. - // - // class TConfidenceLevel - // ---------------------- - // - // Final result of the TLimit algorithm. It is created just after the - // time-consuming part and can be stored in a TFile for further processing. - // It contains light methods to return CLs, CLb and other interesting - // quantities. - // - // The actual algorithm... - // From an input (TLimitDataSource) it produces an output TConfidenceLevel. - // For this, nmc Monte Carlo experiments are performed. - // As usual, the larger this number, the longer the compute time, - // but the better the result. - //Begin_Html - /* - <FONT SIZE=+0> - <p>Supposing that there is a plotfile.root file containing 3 histograms - (signal, background and data), you can imagine doing things like:</p> - <p> - <BLOCKQUOTE><PRE> - TFile* infile=new TFile("plotfile.root","READ"); - infile->cd(); - TH1* sh=(TH1*)infile->Get("signal"); - TH1* bh=(TH1*)infile->Get("background"); - TH1* dh=(TH1*)infile->Get("data"); - TLimitDataSource* mydatasource = new TLimitDataSource(sh,bh,dh); - TConfidenceLevel *myconfidence = TLimit::ComputeLimit(mydatasource,50000); - std::cout << " CLs : " << myconfidence->CLs() << std::endl; - std::cout << " CLsb : " << myconfidence->CLsb() << std::endl; - std::cout << " CLb : " << myconfidence->CLb() << std::endl; - std::cout << "< CLs > : " << myconfidence->GetExpectedCLs_b() << std::endl; - std::cout << "< CLsb > : " << myconfidence->GetExpectedCLsb_b() << std::endl; - std::cout << "< CLb > : " << myconfidence->GetExpectedCLb_b() << std::endl; - delete myconfidence; - delete mydatasource; - infile->Close(); - </PRE></BLOCKQUOTE></p> - <p></p> - <p>More information can still be found on - <a HREF="http://cern.ch/aleph-proj-alphapp/doc/tlimit.html">this</a> page.</p> - </FONT> - */ - //End_Html - // The final object returned... TConfidenceLevel *result = new TConfidenceLevel(nmc); // The random generator used... diff --git a/hist/hist/src/TLimitDataSource.cxx b/hist/hist/src/TLimitDataSource.cxx index e1a84626c2f..42998eb1788 100644 --- a/hist/hist/src/TLimitDataSource.cxx +++ b/hist/hist/src/TLimitDataSource.cxx @@ -2,11 +2,9 @@ // Author: Christophe.Delaere@cern.ch 21/08/2002 /////////////////////////////////////////////////////////////////////////// -// -// TLimitDataSource -// -// This class serves as interface to feed data into the TLimit routines -// +/** \class TLimitDataSource + This class serves as interface to feed data into the TLimit routines +*/ /////////////////////////////////////////////////////////////////////////// #include "TLimitDataSource.h" diff --git a/hist/hist/src/TMultiDimFit.cxx b/hist/hist/src/TMultiDimFit.cxx index cf5f8f85bdf..4cf1377f22e 100644 --- a/hist/hist/src/TMultiDimFit.cxx +++ b/hist/hist/src/TMultiDimFit.cxx @@ -1,41 +1,27 @@ // @(#)root/hist:$Id$ // Author: Christian Holm Christensen 07/11/2000 -//____________________________________________________________________ -// -// -// Begin_Html -/* - </pre> - <H1><A NAME="SECTION00010000000000000000"> - Multidimensional Fits in ROOT</A> - </H1> - - <H1><A NAME="SECTION00020000000000000000"></A> - <A NAME="sec:overview"></A><BR> - Overview - </H1> - - <P> +/** \class TMultiDimFit + \ingroup Hist + + Multidimensional Fits in ROOT. + ## Overview A common problem encountered in different fields of applied science is to find an expression for one physical quantity in terms of several others, which are directly measurable. - <P> An example in high energy physics is the evaluation of the momentum of a charged particle from the observation of its trajectory in a magnetic field. The problem is to relate the momentum of the particle to the observations, which may consists of of positional measurements at intervals along the particle trajectory. - <P> The exact functional relationship between the measured quantities (e.g., the space-points) and the dependent quantity (e.g., the momentum) is in general not known, but one possible way of solving the problem, is to find an expression which reliably approximates the dependence of the momentum on the observations. - <P> This explicit function of the observations can be obtained by a <I>least squares</I> fitting procedure applied to a representive sample of the data, for which the dependent quantity (e.g., momentum) @@ -43,1698 +29,362 @@ used to compute the quantity of interest for new observations of the independent variables. - <P> This class <TT>TMultiDimFit</TT> implements such a procedure in - ROOT. It is largely based on the CERNLIB MUDIFI package - [<A - HREF="TMultiFimFit.html#mudifi">2</A>]. Though the basic concepts are still sound, and + ROOT. It is largely based on the CERNLIB MUDIFI package [2]. + Though the basic concepts are still sound, and therefore kept, a few implementation details have changed, and this - class can take advantage of MINUIT [<A - HREF="TMultiFimFit.html#minuit">4</A>] to improve the errors - of the fitting, thanks to the class <TT>TMinuit</TT>. - - <P> - In [<A - HREF="TMultiFimFit.html#wind72">5</A>] and [<A - HREF="TMultiFimFit.html#wind81">6</A>] H. Wind demonstrates the utility + class can take advantage of MINUIT [4] to improve the errors + of the fitting, thanks to the class TMinuit. + + In [5] and[6] H. Wind demonstrates the utility of this procedure in the context of tracking, magnetic field parameterisation, and so on. The outline of the method used in this class is based on Winds discussion, and I refer these two excellents text for more information. - <P> - And example of usage is given in - <A NAME="tex2html1" - HREF=" - ./examples/multidimfit.C"><TT>$ROOTSYS/tutorials/fit/multidimfit.C</TT></A>. - - <P> - - <H1><A NAME="SECTION00030000000000000000"></A> - <A NAME="sec:method"></A><BR> - The Method - </H1> - - <P> - Let <IMG - WIDTH="18" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img7.gif" - ALT="$ D$"> by the dependent quantity of interest, which depends smoothly - on the observable quantities <!-- MATH - $x_1, \ldots, x_N$ - --> - <IMG - WIDTH="80" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img8.gif" - ALT="$ x_1, \ldots, x_N$">, which we'll denote by - <!-- MATH - $\mathbf{x}$ - --> - <IMG - WIDTH="14" HEIGHT="13" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img9.gif" - ALT="$ \mathbf{x}$">. Given a training sample of <IMG - WIDTH="21" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img10.gif" - ALT="$ M$"> tuples of the form, - (<A NAME="tex2html2" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:AddRow"><TT>TMultiDimFit::AddRow</TT></A>) - <!-- MATH - \begin{displaymath} - \left(\mathbf{x}_j, D_j, E_j\right)\quad, - \end{displaymath} - --> - <P></P><DIV ALIGN="CENTER"> - <IMG - WIDTH="108" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img11.gif" - ALT="$\displaystyle \left(\mathbf{x}_j, D_j, E_j\right)\quad, - $"> - </DIV><P></P> - where <!-- MATH - $\mathbf{x}_j = (x_{1,j},\ldots,x_{N,j})$ - --> - <IMG - WIDTH="148" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img12.gif" - ALT="$ \mathbf{x}_j = (x_{1,j},\ldots,x_{N,j})$"> are <IMG - WIDTH="19" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img13.gif" - ALT="$ N$"> independent - variables, <IMG - WIDTH="24" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img14.gif" - ALT="$ D_j$"> is the known, quantity dependent at <!-- MATH - $\mathbf{x}_j$ - --> - <IMG - WIDTH="20" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img15.gif" - ALT="$ \mathbf{x}_j$">, - and <IMG - WIDTH="23" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img16.gif" - ALT="$ E_j$"> is the square error in <IMG - WIDTH="24" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img14.gif" - ALT="$ D_j$">, the class - <A NAME="tex2html3" - HREF="./TMultiDimFit.html"><TT>TMultiDimFit</TT></A> - will - try to find the parameterization - <P></P> - <DIV ALIGN="CENTER"><A NAME="Dp"></A><!-- MATH - \begin{equation} - D_p(\mathbf{x}) = \sum_{l=1}^{L} c_l \prod_{i=1}^{N} p_{li}\left(x_i\right) - = \sum_{l=1}^{L} c_l F_l(\mathbf{x}) - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="274" HEIGHT="65" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img17.gif" - ALT="$\displaystyle D_p(\mathbf{x}) = \sum_{l=1}^{L} c_l \prod_{i=1}^{N} p_{li}\left(x_i\right) = \sum_{l=1}^{L} c_l F_l(\mathbf{x})$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (1)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> + And example of usage is given in $ROOTSYS/tutorials/fit/multidimfit.C. + + ## The Method + Let \f$ D \f$ by the dependent quantity of interest, which depends smoothly + on the observable quantities \f$ x_1, \ldots, x_N \f$ which we'll denote by + \f$\mathbf{x}\f$. Given a training sample of \f$ M\f$ tuples of the form, (TMultiDimFit::AddRow) + + \f[ + \left(\mathbf{x}_j, D_j, E_j\right)\quad, + \f] + where \f$\mathbf{x}_j = (x_{1,j},\ldots,x_{N,j})\f$ are \f$ N\f$ independent + variables, \f$ D_j\f$ is the known, quantity dependent at \f$\mathbf{x}_j\f$ and \f$ E_j\f$ is + the square error in \f$ D_j\f$, the class will try to find the parameterization + \f[ + D_p(\mathbf{x}) = \sum_{l=1}^{L} c_l \prod_{i=1}^{N} p_{li}\left(x_i\right) + = \sum_{l=1}^{L} c_l F_l(\mathbf{x}) + \f] such that - <P></P> - <DIV ALIGN="CENTER"><A NAME="S"></A><!-- MATH - \begin{equation} - S \equiv \sum_{j=1}^{M} \left(D_j - D_p\left(\mathbf{x}_j\right)\right)^2 - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="172" HEIGHT="65" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img18.gif" - ALT="$\displaystyle S \equiv \sum_{j=1}^{M} \left(D_j - D_p\left(\mathbf{x}_j\right)\right)^2$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (2)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - is minimal. Here <!-- MATH - $p_{li}(x_i)$ - --> - <IMG - WIDTH="48" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img19.gif" - ALT="$ p_{li}(x_i)$"> are monomials, or Chebyshev or Legendre - polynomials, labelled <!-- MATH - $l = 1, \ldots, L$ - --> - <IMG - WIDTH="87" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img20.gif" - ALT="$ l = 1, \ldots, L$">, in each variable - <IMG - WIDTH="18" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img21.gif" - ALT="$ x_i$">, <!-- MATH - $i=1, \ldots, N$ - --> - <IMG - WIDTH="91" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img22.gif" - ALT="$ i=1, \ldots, N$">. - - <P> - So what <TT>TMultiDimFit</TT> does, is to determine the number of - terms <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img23.gif" - ALT="$ L$">, and then <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img23.gif" - ALT="$ L$"> terms (or functions) <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$">, and the <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img23.gif" - ALT="$ L$"> - coefficients <IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img25.gif" - ALT="$ c_l$">, so that <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> is minimal - (<A NAME="tex2html4" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:FindParameterization"><TT>TMultiDimFit::FindParameterization</TT></A>). - - <P> - Of course it's more than a little unlikely that <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> will ever become + + \f[ + S \equiv \sum_{j=1}^{M} \left(D_j - D_p\left(\mathbf{x}_j\right)\right)^2 + \f] + is minimal. Here \f$p_{li}(x_i)\f$ are monomials, or Chebyshev or Legendre + polynomials, labelled \f$l = 1, \ldots, L\f$, in each variable \f$ x_i\f$,\f$ i=1, \ldots, N\f$. + + So what TMultiDimFit does, is to determine the number of terms \f$ L\f$, and then \f$ L\f$ terms + (or functions) \f$ F_l\f$, and the \f$ L\f$ coefficients \f$ c_l\f$, so that \f$ S\f$ is minimal + (TMultiDimFit::FindParameterization). + + Of course it's more than a little unlikely that \f$ S\f$ will ever become exact zero as a result of the procedure outlined below. Therefore, the - user is asked to provide a minimum relative error <IMG - WIDTH="11" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img27.gif" - ALT="$ \epsilon$"> - (<A NAME="tex2html5" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:SetMinRelativeError"><TT>TMultiDimFit::SetMinRelativeError</TT></A>), and <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> - will be considered minimized when - <!-- MATH - \begin{displaymath} - R = \frac{S}{\sum_{j=1}^M D_j^2} < \epsilon - \end{displaymath} - --> - <P></P><DIV ALIGN="CENTER"> - <IMG - WIDTH="132" HEIGHT="51" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img28.gif" - ALT="$\displaystyle R = \frac{S}{\sum_{j=1}^M D_j^2} < \epsilon - $"> - </DIV><P></P> - - <P> + user is asked to provide a minimum relative error \f$ \epsilon\f$ (TMultiDimFit::SetMinRelativeError), + and \f$ S\f$ will be considered minimized when + + \f[ + R = \frac{S}{\sum_{j=1}^M D_j^2} < \epsilon + \f] Optionally, the user may impose a functional expression by specifying - the powers of each variable in <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img23.gif" - ALT="$ L$"> specified functions <!-- MATH - $F_1, \ldots, - F_L$ - --> - <IMG - WIDTH="79" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img29.gif" - ALT="$ F_1, \ldots, - F_L$"> (<A NAME="tex2html6" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:SetPowers"><TT>TMultiDimFit::SetPowers</TT></A>). In that case, only the - coefficients <IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img25.gif" - ALT="$ c_l$"> is calculated by the class. - - <P> - - <H2><A NAME="SECTION00031000000000000000"></A> - <A NAME="sec:selection"></A><BR> - Limiting the Number of Terms - </H2> - - <P> - As always when dealing with fits, there's a real chance of - <I>over fitting</I>. As is well-known, it's always possible to fit an - <IMG - WIDTH="46" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img30.gif" - ALT="$ N-1$"> polynomial in <IMG - WIDTH="13" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img31.gif" - ALT="$ x$"> to <IMG - WIDTH="19" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img13.gif" - ALT="$ N$"> points <IMG - WIDTH="41" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img32.gif" - ALT="$ (x,y)$"> with <!-- MATH - $\chi^2 = 0$ - --> - <IMG - WIDTH="50" HEIGHT="33" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img33.gif" - ALT="$ \chi^2 = 0$">, but - the polynomial is not likely to fit new data at all - [<A - HREF="TMultiFimFit.html#bevington">1</A>]. Therefore, the user is asked to provide an upper - limit, <IMG - WIDTH="41" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img34.gif" - ALT="$ L_{max}$"> to the number of terms in <IMG - WIDTH="25" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img35.gif" - ALT="$ D_p$"> - (<A NAME="tex2html7" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:SetMaxTerms"><TT>TMultiDimFit::SetMaxTerms</TT></A>). - - <P> - However, since there's an infinite number of <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$"> to choose from, the - user is asked to give the maximum power. <IMG - WIDTH="49" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img36.gif" - ALT="$ P_{max,i}$">, of each variable - <IMG - WIDTH="18" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img21.gif" - ALT="$ x_i$"> to be considered in the minimization of <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> - (<A NAME="tex2html8" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:SetMaxPowers"><TT>TMultiDimFit::SetMaxPowers</TT></A>). - - <P> - One way of obtaining values for the maximum power in variable <IMG - WIDTH="10" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img37.gif" - ALT="$ i$">, is - to perform a regular fit to the dependent quantity <IMG - WIDTH="18" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img7.gif" - ALT="$ D$">, using a - polynomial only in <IMG - WIDTH="18" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img21.gif" - ALT="$ x_i$">. The maximum power is <IMG - WIDTH="49" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img36.gif" - ALT="$ P_{max,i}$"> is then the + the powers of each variable in \f$ L\f$ specified functions \f$ F_1, \ldots,F_L\f$ (TMultiDimFit::SetPowers). + In that case, only the coefficients \f$ c_l\f$ is calculated by the class. + + ## Limiting the Number of Terms + As always when dealing with fits, there's a real chance of *over fitting*. As is well-known, it's + always possible to fit an \f$ N-1\f$ polynomial in \f$ x\f$ to \f$ N\f$ points \f$ (x,y)\f$ with + \f$\chi^2 = 0\f$, but the polynomial is not likely to fit new data at all [1]. + Therefore, the user is asked to provide an upper limit, \f$ L_{max}\f$ to the number of terms in + \f$ D_p\f$ (TMultiDimFit::SetMaxTerms). + + However, since there's an infinite number of \f$ F_l\f$ to choose from, the + user is asked to give the maximum power. \f$ P_{max,i}\f$, of each variable + \f$ x_i\f$ to be considered in the minimization of \f$ S\f$ (TMultiDimFit::SetMaxPowers). + + One way of obtaining values for the maximum power in variable \f$ i\f$, is + to perform a regular fit to the dependent quantity \f$ D\f$, using a + polynomial only in \f$ x_i\f$. The maximum power is \f$ P_{max,i}\f$ is then the power that does not significantly improve the one-dimensional - least-square fit over <IMG - WIDTH="18" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img21.gif" - ALT="$ x_i$"> to <IMG - WIDTH="18" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img7.gif" - ALT="$ D$"> [<A - HREF="TMultiFimFit.html#wind72">5</A>]. - - <P> - There are still a huge amount of possible choices for <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$">; in fact - there are <!-- MATH - $\prod_{i=1}^{N} (P_{max,i} + 1)$ - --> - <IMG - WIDTH="125" HEIGHT="39" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img38.gif" - ALT="$ \prod_{i=1}^{N} (P_{max,i} + 1)$"> possible + least-square fit over \f$ x_i\f$ to \f$ D\f$ [5]. + + There are still a huge amount of possible choices for \f$ F_l\f$; in fact + there are \f$\prod_{i=1}^{N} (P_{max,i} + 1)\f$ possible choices. Obviously we need to limit this. To this end, the user is - asked to set a <I>power control limit</I>, <IMG - WIDTH="17" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img39.gif" - ALT="$ Q$"> - (<A NAME="tex2html9" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:SetPowerLimit"><TT>TMultiDimFit::SetPowerLimit</TT></A>), and a function - <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$"> is only accepted if - <!-- MATH - \begin{displaymath} - Q_l = \sum_{i=1}^{N} \frac{P_{li}}{P_{max,i}} < Q - \end{displaymath} - --> - <P></P><DIV ALIGN="CENTER"> - <IMG - WIDTH="151" HEIGHT="65" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img40.gif" - ALT="$\displaystyle Q_l = \sum_{i=1}^{N} \frac{P_{li}}{P_{max,i}} < Q - $"> - </DIV><P></P> - where <IMG - WIDTH="24" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img41.gif" - ALT="$ P_{li}$"> is the leading power of variable <IMG - WIDTH="18" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img21.gif" - ALT="$ x_i$"> in function - <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$">. (<A NAME="tex2html10" - HREF=" - - ./TMultiDimFit.html#TMultiDimFit:MakeCandidates"><TT>TMultiDimFit::MakeCandidates</TT></A>). So the number of - functions increase with <IMG - WIDTH="17" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img39.gif" - ALT="$ Q$"> (1, 2 is fine, 5 is way out). - - <P> - - <H2><A NAME="SECTION00032000000000000000"> - Gram-Schmidt Orthogonalisation</A> - </H2> - - <P> + asked to set a *power control limit*, \f$ Q\f$ (TMultiDimFit::SetPowerLimit), and a function + \f$ F_l\f$ is only accepted if + \f[ + Q_l = \sum_{i=1}^{N} \frac{P_{li}}{P_{max,i}} < Q + \f] + where \f$ P_{li}\f$ is the leading power of variable \f$ x_i\f$ in function \f$ F_l\f$ (TMultiDimFit::MakeCandidates). + So the number of functions increase with \f$ Q\f$ (1, 2 is fine, 5 is way out). + + ## Gram-Schmidt Orthogonalisation</A> To further reduce the number of functions in the final expression, - only those functions that significantly reduce <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> is chosen. What + only those functions that significantly reduce \f$ S\f$ is chosen. What `significant' means, is chosen by the user, and will be - discussed below (see <A HREF="TMultiFimFit.html#sec:selectiondetail">2.3</A>). - - <P> - The functions <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$"> are generally not orthogonal, which means one will - have to evaluate all possible <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$">'s over all data-points before - finding the most significant [<A - HREF="TMultiFimFit.html#bevington">1</A>]. We can, however, do - better then that. By applying the <I>modified Gram-Schmidt - orthogonalisation</I> algorithm [<A - HREF="TMultiFimFit.html#wind72">5</A>] [<A - HREF="TMultiFimFit.html#golub">3</A>] to the - functions <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$">, we can evaluate the contribution to the reduction of - <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> from each function in turn, and we may delay the actual inversion - of the curvature-matrix - (<A NAME="tex2html11" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:MakeGramSchmidt"><TT>TMultiDimFit::MakeGramSchmidt</TT></A>). - - <P> - So we are let to consider an <IMG - WIDTH="52" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img42.gif" - ALT="$ M\times L$"> matrix <!-- MATH - $\mathsf{F}$ - --> - <IMG - WIDTH="13" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img43.gif" - ALT="$ \mathsf{F}$">, an + discussed below (see [2.3](TMultiFimFit.html#sec:selectiondetail)). + + The functions \f$ F_l\f$ are generally not orthogonal, which means one will + have to evaluate all possible \f$ F_l\f$'s over all data-points before + finding the most significant [1]. We can, however, do + better then that. By applying the *modified Gram-Schmidt + orthogonalisation* algorithm [5] [3] to the + functions \f$ F_l\f$, we can evaluate the contribution to the reduction of + \f$ S\f$ from each function in turn, and we may delay the actual inversion + of the curvature-matrix (TMultiDimFit::MakeGramSchmidt). + + So we are let to consider an \f$ M\times L\f$ matrix \f$\mathsf{F}\f$, an element of which is given by - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:Felem"></A><!-- MATH - \begin{equation} - f_{jl} = F_j\left(x_{1j} , x_{2j}, \ldots, x_{Nj}\right) - = F_l(\mathbf{x}_j)\, \quad\mbox{with}~j=1,2,\ldots,M, - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="260" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img44.gif" - ALT="$\displaystyle f_{jl} = F_j\left(x_{1j} , x_{2j}, \ldots, x_{Nj}\right) = F_l(\mathbf{x}_j) $"> with<IMG - WIDTH="120" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img45.gif" - ALT="$\displaystyle j=1,2,\ldots,M,$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (3)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - where <IMG - WIDTH="12" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img46.gif" - ALT="$ j$"> labels the <IMG - WIDTH="21" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img10.gif" - ALT="$ M$"> rows in the training sample and <IMG - WIDTH="9" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img47.gif" - ALT="$ l$"> labels - <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img23.gif" - ALT="$ L$"> functions of <IMG - WIDTH="19" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img13.gif" - ALT="$ N$"> variables, and <IMG - WIDTH="53" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img48.gif" - ALT="$ L \leq M$">. That is, <IMG - WIDTH="23" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img49.gif" - ALT="$ f_{jl}$"> is - the term (or function) numbered <IMG - WIDTH="9" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img47.gif" - ALT="$ l$"> evaluated at the data point - <IMG - WIDTH="12" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img46.gif" - ALT="$ j$">. We have to normalise <!-- MATH - $\mathbf{x}_j$ - --> - <IMG - WIDTH="20" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img15.gif" - ALT="$ \mathbf{x}_j$"> to <IMG - WIDTH="48" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img50.gif" - ALT="$ [-1,1]$"> for this to - succeed [<A - HREF="TMultiFimFit.html#wind72">5</A>] - (<A NAME="tex2html12" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:MakeNormalized"><TT>TMultiDimFit::MakeNormalized</TT></A>). We then define a - matrix <!-- MATH - $\mathsf{W}$ - --> - <IMG - WIDTH="19" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img51.gif" - ALT="$ \mathsf{W}$"> of which the columns <!-- MATH - $\mathbf{w}_j$ - --> - <IMG - WIDTH="24" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img52.gif" - ALT="$ \mathbf{w}_j$"> are given by - <BR> - <DIV ALIGN="CENTER"><A NAME="eq:wj"></A><!-- MATH - \begin{eqnarray} - \mathbf{w}_1 &=& \mathbf{f}_1 = F_1\left(\mathbf x_1\right)\\ - \mathbf{w}_l &=& \mathbf{f}_l - \sum^{l-1}_{k=1} \frac{\mathbf{f}_l \bullet - \mathbf{w}_k}{\mathbf{w}_k^2}\mathbf{w}_k\,. - \end{eqnarray} - --> - <TABLE CELLPADDING="0" ALIGN="CENTER" WIDTH="100%"> - <TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"><IMG - WIDTH="25" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img53.gif" - ALT="$\displaystyle \mathbf{w}_1$"></TD> - <TD WIDTH="10" ALIGN="CENTER" NOWRAP><IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img54.gif" - ALT="$\displaystyle =$"></TD> - <TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="87" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img55.gif" - ALT="$\displaystyle \mathbf{f}_1 = F_1\left(\mathbf x_1\right)$"></TD> - <TD WIDTH=10 ALIGN="RIGHT"> - (4)</TD></TR> - <TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"><IMG - WIDTH="22" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img56.gif" - ALT="$\displaystyle \mathbf{w}_l$"></TD> - <TD WIDTH="10" ALIGN="CENTER" NOWRAP><IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img54.gif" - ALT="$\displaystyle =$"></TD> - <TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="138" HEIGHT="66" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img57.gif" - ALT="$\displaystyle \mathbf{f}_l - \sum^{l-1}_{k=1} \frac{\mathbf{f}_l \bullet - \mathbf{w}_k}{\mathbf{w}_k^2}\mathbf{w}_k .$"></TD> - <TD WIDTH=10 ALIGN="RIGHT"> - (5)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - and <!-- MATH - $\mathbf{w}_{l}$ - --> - <IMG - WIDTH="22" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img58.gif" - ALT="$ \mathbf{w}_{l}$"> is the component of <!-- MATH - $\mathbf{f}_{l}$ - --> - <IMG - WIDTH="15" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img59.gif" - ALT="$ \mathbf{f}_{l}$"> orthogonal - to <!-- MATH - $\mathbf{w}_{1}, \ldots, \mathbf{w}_{l-1}$ - --> - <IMG - WIDTH="97" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img60.gif" - ALT="$ \mathbf{w}_{1}, \ldots, \mathbf{w}_{l-1}$">. Hence we obtain - [<A - HREF="TMultiFimFit.html#golub">3</A>], - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:worto"></A><!-- MATH - \begin{equation} - \mathbf{w}_k\bullet\mathbf{w}_l = 0\quad\mbox{if}~k \neq l\quad. - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="87" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img61.gif" - ALT="$\displaystyle \mathbf{w}_k\bullet\mathbf{w}_l = 0$"> if<IMG - WIDTH="65" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img62.gif" - ALT="$\displaystyle k \neq l\quad.$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (6)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - - <P> - We now take as a new model <!-- MATH - $\mathsf{W}\mathbf{a}$ - --> - <IMG - WIDTH="28" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img63.gif" - ALT="$ \mathsf{W}\mathbf{a}$">. We thus want to + \f[ + f_{jl} = F_j\left(x_{1j} , x_{2j}, \ldots, x_{Nj}\right) + = F_l(\mathbf{x}_j)\, \quad\mbox{with}~j=1,2,\ldots,M, + \f] + where \f$ j\f$ labels the \f$ M\f$ rows in the training sample and \f$ l\f$ labels + \f$ L\f$ functions of \f$ N\f$ variables, and \f$ L \leq M\f$. That is, \f$ f_{jl}\f$ is + the term (or function) numbered \f$ l\f$ evaluated at the data point + \f$ j\f$. We have to normalise \f$\mathbf{x}_j\f$ to \f$ [-1,1]\f$ for this to + succeed [5] (TMultiDimFit::MakeNormalized). We then define a + matrix \f$\mathsf{W}\f$ of which the columns \f$\mathbf{w}_j\f$ are given by + \f{eqnarray*}{ + \mathbf{w}_1 &=& \mathbf{f}_1 = F_1\left(\mathbf x_1\right)\\ + \mathbf{w}_l &=& \mathbf{f}_l - \sum^{l-1}_{k=1} \frac{\mathbf{f}_l \bullet + \mathbf{w}_k}{\mathbf{w}_k^2}\mathbf{w}_k\,. + \f} + and \f$\mathbf{w}_{l}\f$ is the component of \f$\mathbf{f}_{l} \f$ orthogonal + to \f$\mathbf{w}_{1}, \ldots, \mathbf{w}_{l-1}\f$. Hence we obtain [3], + \f[ + \mathbf{w}_k\bullet\mathbf{w}_l = 0\quad\mbox{if}~k \neq l\quad. + \f] + We now take as a new model \f$\mathsf{W}\mathbf{a}\f$. We thus want to minimize - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:S"></A><!-- MATH - \begin{equation} - S\equiv \left(\mathbf{D} - \mathsf{W}\mathbf{a}\right)^2\quad, - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="136" HEIGHT="38" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img64.gif" - ALT="$\displaystyle S\equiv \left(\mathbf{D} - \mathsf{W}\mathbf{a}\right)^2\quad,$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (7)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - where <!-- MATH - $\mathbf{D} = \left(D_1,\ldots,D_M\right)$ - --> - <IMG - WIDTH="137" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img65.gif" - ALT="$ \mathbf{D} = \left(D_1,\ldots,D_M\right)$"> is a vector of the + \f[ + S\equiv \left(\mathbf{D} - \mathsf{W}\mathbf{a}\right)^2\quad, + \f] + where \f$\mathbf{D} = \left(D_1,\ldots,D_M\right)\f$ is a vector of the dependent quantity in the sample. Differentiation with respect to - <IMG - WIDTH="19" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img66.gif" - ALT="$ a_j$"> gives, using (<A HREF="TMultiFimFit.html#eq:worto">6</A>), - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:dS"></A><!-- MATH - \begin{equation} - \mathbf{D}\bullet\mathbf{w}_l - a_l\mathbf{w}_l^2 = 0 - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="134" HEIGHT="35" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img67.gif" - ALT="$\displaystyle \mathbf{D}\bullet\mathbf{w}_l - a_l\mathbf{w}_l^2 = 0$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (8)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> + \f$ a_j\f$ gives, using [6], + \f[ + \mathbf{D}\bullet\mathbf{w}_l - a_l\mathbf{w}_l^2 = 0 + \f] or - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:dS2"></A><!-- MATH - \begin{equation} - a_l = \frac{\mathbf{D}_l\bullet\mathbf{w}_l}{\mathbf{w}_l^2} - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="95" HEIGHT="51" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img68.gif" - ALT="$\displaystyle a_l = \frac{\mathbf{D}_l\bullet\mathbf{w}_l}{\mathbf{w}_l^2}$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (9)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - Let <IMG - WIDTH="21" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img69.gif" - ALT="$ S_j$"> be the sum of squares of residuals when taking <IMG - WIDTH="12" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img46.gif" - ALT="$ j$"> functions + \f[ + a_l = \frac{\mathbf{D}_l\bullet\mathbf{w}_l}{\mathbf{w}_l^2} + \f] + Let \f$ S_j\f$ be the sum of squares of residuals when taking \f$ j\f$ functions into account. Then - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:Sj"></A><!-- MATH - \begin{equation} - S_l = \left[\mathbf{D} - \sum^l_{k=1} a_k\mathbf{w}_k\right]^2 - = \mathbf{D}^2 - 2\mathbf{D} \sum^l_{k=1} a_k\mathbf{w}_k - + \sum^l_{k=1} a_k^2\mathbf{w}_k^2 - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="394" HEIGHT="72" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img70.gif" - ALT="$\displaystyle S_l = \left[\mathbf{D} - \sum^l_{k=1} a_k\mathbf{w}_k\right]^2 = ... - ...2 - 2\mathbf{D} \sum^l_{k=1} a_k\mathbf{w}_k + \sum^l_{k=1} a_k^2\mathbf{w}_k^2$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (10)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - Using (<A HREF="TMultiFimFit.html#eq:dS2">9</A>), we see that - <BR> - <DIV ALIGN="CENTER"><A NAME="eq:sj2"></A><!-- MATH - \begin{eqnarray} - S_l &=& \mathbf{D}^2 - 2 \sum^l_{k=1} a_k^2\mathbf{w}_k^2 + - \sum^j_{k=1} a_k^2\mathbf{w}_k^2\nonumber\\ - &=& \mathbf{D}^2 - \sum^l_{k=1} a_k^2\mathbf{w}_k^2\nonumber\\ - &=& \mathbf{D}^2 - \sum^l_{k=1} \frac{\left(\mathbf D\bullet \mathbf - w_k\right)}{\mathbf w_k^2} - \end{eqnarray} - --> - <TABLE CELLPADDING="0" ALIGN="CENTER" WIDTH="100%"> - <TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"><IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img71.gif" - ALT="$\displaystyle S_l$"></TD> - <TD WIDTH="10" ALIGN="CENTER" NOWRAP><IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img54.gif" - ALT="$\displaystyle =$"></TD> - <TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="201" HEIGHT="67" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img72.gif" - ALT="$\displaystyle \mathbf{D}^2 - 2 \sum^l_{k=1} a_k^2\mathbf{w}_k^2 + - \sum^j_{k=1} a_k^2\mathbf{w}_k^2$"></TD> - <TD WIDTH=10 ALIGN="RIGHT"> - </TD></TR> - <TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"> </TD> - <TD WIDTH="10" ALIGN="CENTER" NOWRAP><IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img54.gif" - ALT="$\displaystyle =$"></TD> - <TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="108" HEIGHT="66" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img73.gif" - ALT="$\displaystyle \mathbf{D}^2 - \sum^l_{k=1} a_k^2\mathbf{w}_k^2$"></TD> - <TD WIDTH=10 ALIGN="RIGHT"> - </TD></TR> - <TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"> </TD> - <TD WIDTH="10" ALIGN="CENTER" NOWRAP><IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img54.gif" - ALT="$\displaystyle =$"></TD> - <TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="137" HEIGHT="66" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img74.gif" - ALT="$\displaystyle \mathbf{D}^2 - \sum^l_{k=1} \frac{\left(\mathbf D\bullet \mathbf - w_k\right)}{\mathbf w_k^2}$"></TD> - <TD WIDTH=10 ALIGN="RIGHT"> - (11)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - - <P> - So for each new function <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$"> included in the model, we get a - reduction of the sum of squares of residuals of <!-- MATH - $a_l^2\mathbf{w}_l^2$ - --> - <IMG - WIDTH="40" HEIGHT="33" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img75.gif" - ALT="$ a_l^2\mathbf{w}_l^2$">, - where <!-- MATH - $\mathbf{w}_l$ - --> - <IMG - WIDTH="22" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img76.gif" - ALT="$ \mathbf{w}_l$"> is given by (<A HREF="TMultiFimFit.html#eq:wj">4</A>) and <IMG - WIDTH="17" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img77.gif" - ALT="$ a_l$"> by - (<A HREF="TMultiFimFit.html#eq:dS2">9</A>). Thus, using the Gram-Schmidt orthogonalisation, we + \f[ + S_l = \left[\mathbf{D} - \sum^l_{k=1} a_k\mathbf{w}_k\right]^2 + = \mathbf{D}^2 - 2\mathbf{D} \sum^l_{k=1} a_k\mathbf{w}_k + + \sum^l_{k=1} a_k^2\mathbf{w}_k^2 + \f] + Using [9], we see that + \f[ + S_l &=& \mathbf{D}^2 - 2 \sum^l_{k=1} a_k^2\mathbf{w}_k^2 + + \sum^j_{k=1} a_k^2\mathbf{w}_k^2\nonumber\\ + &=& \mathbf{D}^2 - \sum^l_{k=1} a_k^2\mathbf{w}_k^2\nonumber\\ + &=& \mathbf{D}^2 - \sum^l_{k=1} \frac{\left(\mathbf D\bullet \mathbf + w_k\right)}{\mathbf w_k^2} + \f] + So for each new function \f$ F_l\f$ included in the model, we get a + reduction of the sum of squares of residuals of \f$a_l^2\mathbf{w}_l^2\f$, + where \f$\mathbf{w}_l\f$ is given by [4] and \f$ a_l\f$ by [9]. Thus, using + the Gram-Schmidt orthogonalisation, we can decide if we want to include this function in the final model, - <I>before</I> the matrix inversion. - - <P> - - <H2><A NAME="SECTION00033000000000000000"></A> - <A NAME="sec:selectiondetail"></A><BR> - Function Selection Based on Residual - </H2> - - <P> - Supposing that <IMG - WIDTH="42" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img78.gif" - ALT="$ L-1$"> steps of the procedure have been performed, the - problem now is to consider the <!-- MATH - $L^{\mbox{th}}$ - --> - <IMG - WIDTH="31" HEIGHT="20" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img79.gif" - ALT="$ L^{\mbox{th}}$"> function. - - <P> + *before* the matrix inversion. + + ## Function Selection Based on Residual + Supposing that \f$ L-1\f$ steps of the procedure have been performed, the + problem now is to consider the \f$L^{\mbox{th}}\f$ function. + The sum of squares of residuals can be written as - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:sums"></A><!-- MATH - \begin{equation} - S_L = \textbf{D}^T\bullet\textbf{D} - - \sum^L_{l=1}a^2_l\left(\textbf{w}_l^T\bullet\textbf{w}_l\right) - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="232" HEIGHT="65" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img80.gif" - ALT="$\displaystyle S_L = \textbf{D}^T\bullet\textbf{D} - \sum^L_{l=1}a^2_l\left(\textbf{w}_l^T\bullet\textbf{w}_l\right)$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (12)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - where the relation (<A HREF="TMultiFimFit.html#eq:dS2">9</A>) have been taken into account. The - contribution of the <!-- MATH - $L^{\mbox{th}}$ - --> - <IMG - WIDTH="31" HEIGHT="20" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img79.gif" - ALT="$ L^{\mbox{th}}$"> function to the reduction of S, is + \f[ + S_L = \textbf{D}^T\bullet\textbf{D} - + \sum^L_{l=1}a^2_l\left(\textbf{w}_l^T\bullet\textbf{w}_l\right) + \f] + where the relation [9] have been taken into account. The + contribution of the \f$L^{\mbox{th}}\f$ function to the reduction of S, is given by - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:dSN"></A><!-- MATH - \begin{equation} - \Delta S_L = a^2_L\left(\textbf{w}_L^T\bullet\textbf{w}_L\right) - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="154" HEIGHT="36" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img81.gif" - ALT="$\displaystyle \Delta S_L = a^2_L\left(\textbf{w}_L^T\bullet\textbf{w}_L\right)$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (13)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - - <P> - Two test are now applied to decide whether this <!-- MATH - $L^{\mbox{th}}$ - --> - <IMG - WIDTH="31" HEIGHT="20" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img79.gif" - ALT="$ L^{\mbox{th}}$"> + \f[ + \Delta S_L = a^2_L\left(\textbf{w}_L^T\bullet\textbf{w}_L\right) + \f] + Two test are now applied to decide whether this \f$L^{\mbox{th}}\f$ function is to be included in the final expression, or not. - <P> - - <H3><A NAME="SECTION00033100000000000000"></A> - <A NAME="testone"></A><BR> - Test 1 - </H3> - - <P> - Denoting by <IMG - WIDTH="43" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img82.gif" - ALT="$ H_{L-1}$"> the subspace spanned by - <!-- MATH - $\textbf{w}_1,\ldots,\textbf{w}_{L-1}$ - --> - <IMG - WIDTH="102" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img83.gif" - ALT="$ \textbf{w}_1,\ldots,\textbf{w}_{L-1}$"> the function <!-- MATH - $\textbf{w}_L$ - --> - <IMG - WIDTH="27" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img5.gif" - ALT="$ \textbf {w}_L$"> is - by construction (see (<A HREF="TMultiFimFit.html#eq:wj">4</A>)) the projection of the function - <IMG - WIDTH="24" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img84.gif" - ALT="$ F_L$"> onto the direction perpendicular to <IMG - WIDTH="43" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img82.gif" - ALT="$ H_{L-1}$">. Now, if the - length of <!-- MATH - $\textbf{w}_L$ - --> - <IMG - WIDTH="27" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img5.gif" - ALT="$ \textbf {w}_L$"> (given by <!-- MATH - $\textbf{w}_L\bullet\textbf{w}_L$ - --> - <IMG - WIDTH="65" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img85.gif" - ALT="$ \textbf{w}_L\bullet\textbf{w}_L$">) - is very small compared to the length of <!-- MATH - $\textbf{f}_L$ - --> - <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img3.gif" - ALT="$ \textbf {f}_L$"> this new + ## Test 1 + Denoting by \f$ H_{L-1}\f$ the subspace spanned by \f$\textbf{w}_1,\ldots,\textbf{w}_{L-1}\f$ + the function \d$\textbf{w}_L\d$ is by construction (see 4) the projection of the function + \f$ F_L\f$ onto the direction perpendicular to \f$ H_{L-1}\f$. Now, if the + length of \f$\textbf{w}_L\f$ (given by \f$\textbf{w}_L\bullet\textbf{w}_L\f$) + is very small compared to the length of \f$\textbf{f}_L\f$ this new function can not contribute much to the reduction of the sum of squares of residuals. The test consists then in calculating the angle - <IMG - WIDTH="12" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img1.gif" - ALT="$ \theta $"> between the two vectors <!-- MATH - $\textbf{w}_L$ - --> - <IMG - WIDTH="27" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img5.gif" - ALT="$ \textbf {w}_L$"> and <!-- MATH - $\textbf{f}_L$ - --> - <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img3.gif" - ALT="$ \textbf {f}_L$"> - (see also figure <A HREF="TMultiFimFit.html#fig:thetaphi">1</A>) and requiring that it's - <I>greater</I> then a threshold value which the user must set - (<A NAME="tex2html14" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:SetMinAngle"><TT>TMultiDimFit::SetMinAngle</TT></A>). - - <P> - - <P></P> - <DIV ALIGN="CENTER"><A NAME="fig:thetaphi"></A><A NAME="519"></A> - <TABLE> - <CAPTION ALIGN="BOTTOM"><STRONG>Figure 1:</STRONG> - (a) Angle <IMG - WIDTH="12" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img1.gif" - ALT="$ \theta $"> between <!-- MATH - $\textbf{w}_l$ - --> - <IMG - WIDTH="22" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img2.gif" - ALT="$ \textbf {w}_l$"> and - <!-- MATH - $\textbf{f}_L$ - --> - <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img3.gif" - ALT="$ \textbf {f}_L$">, (b) angle <IMG - WIDTH="14" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img4.gif" - ALT="$ \phi $"> between <!-- MATH - $\textbf{w}_L$ - --> - <IMG - WIDTH="27" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img5.gif" - ALT="$ \textbf {w}_L$"> and - <!-- MATH - $\textbf{D}$ - --> - <IMG - WIDTH="18" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img6.gif" - ALT="$ \textbf {D}$"></CAPTION> - <TR><TD><IMG - WIDTH="466" HEIGHT="172" BORDER="0" - SRC="gif/multidimfit_img86.gif" - ALT="\begin{figure}\begin{center} - \begin{tabular}{p{.4\textwidth}p{.4\textwidth}} - \... - ... \put(80,100){$\mathbf{D}$} - \end{picture} \end{tabular} \end{center}\end{figure}"></TD></TR> - </TABLE> - </DIV><P></P> - - <P> - - <H3><A NAME="SECTION00033200000000000000"></A> <A NAME="testtwo"></A><BR> - Test 2 - </H3> - - <P> - Let <!-- MATH - $\textbf{D}$ - --> - <IMG - WIDTH="18" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img6.gif" - ALT="$ \textbf {D}$"> be the data vector to be fitted. As illustrated in - figure <A HREF="TMultiFimFit.html#fig:thetaphi">1</A>, the <!-- MATH - $L^{\mbox{th}}$ - --> - <IMG - WIDTH="31" HEIGHT="20" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img79.gif" - ALT="$ L^{\mbox{th}}$"> function <!-- MATH - $\textbf{w}_L$ - --> - <IMG - WIDTH="27" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img5.gif" - ALT="$ \textbf {w}_L$"> - will contribute significantly to the reduction of <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$">, if the angle - <!-- MATH - $\phi^\prime$ - --> - <IMG - WIDTH="18" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img87.gif" - ALT="$ \phi^\prime$"> between <!-- MATH - $\textbf{w}_L$ - --> - <IMG - WIDTH="27" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img5.gif" - ALT="$ \textbf {w}_L$"> and <!-- MATH - $\textbf{D}$ - --> - <IMG - WIDTH="18" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img6.gif" - ALT="$ \textbf {D}$"> is smaller than - an upper limit <IMG - WIDTH="14" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img4.gif" - ALT="$ \phi $">, defined by the user - (<A NAME="tex2html15" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:SetMaxAngle"><TT>TMultiDimFit::SetMaxAngle</TT></A>) - - <P> + \f$ \theta \f$ between the two vectors \f$\textbf{w}_L\f$ \f$ \textbf {f}_L\f$ + (see also figure 1) and requiring that it's + *greater* then a threshold value which the user must set (TMultiDimFit::SetMinAngle). + + \image html multidimfit_img86.gif "Figure 1: (a) angle \f$\theta\f$ between \f$\textbf{w}_l\f$ and \f$\textbf{f}_L\f$, (b) angle \f$ \phi \f$ between \f$\textbf{w}_L\f$ and \f$\textbf{D}\f$" + + ## Test 2 + Let \f$\textbf{D}\f$ be the data vector to be fitted. As illustrated in + figure 1, the \f$L^{\mbox{th}}\f$ function \f$\textbf{w}_L\f$ + will contribute significantly to the reduction of \f$ S\f$, if the angle + \f$\phi^\prime\f$ between \f$\textbf{w}_L\f$ and \f$\textbf{D}\f$ is smaller than + an upper limit \f$ \phi \f$, defined by the user (MultiDimFit::SetMaxAngle) + However, the method automatically readjusts the value of this angle while fitting is in progress, in order to make the selection criteria less and less difficult to be fulfilled. The result is that the - functions contributing most to the reduction of <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> are chosen first - (<A NAME="tex2html16" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:TestFunction"><TT>TMultiDimFit::TestFunction</TT></A>). - - <P> - In case <IMG - WIDTH="14" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img4.gif" - ALT="$ \phi $"> isn't defined, an alternative method of - performing this second test is used: The <!-- MATH - $L^{\mbox{th}}$ - --> - <IMG - WIDTH="31" HEIGHT="20" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img79.gif" - ALT="$ L^{\mbox{th}}$"> function - <!-- MATH - $\textbf{f}_L$ - --> - <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img3.gif" - ALT="$ \textbf {f}_L$"> is accepted if (refer also to equation (<A HREF="TMultiFimFit.html#eq:dSN">13</A>)) - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:dSN2"></A><!-- MATH - \begin{equation} - \Delta S_L > \frac{S_{L-1}}{L_{max}-L} - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="129" HEIGHT="51" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img88.gif" - ALT="$\displaystyle \Delta S_L > \frac{S_{L-1}}{L_{max}-L}$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (14)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - where <IMG - WIDTH="40" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img89.gif" - ALT="$ S_{L-1}$"> is the sum of the <IMG - WIDTH="42" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img78.gif" - ALT="$ L-1$"> first residuals from the - <IMG - WIDTH="42" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img78.gif" - ALT="$ L-1$"> functions previously accepted; and <IMG - WIDTH="41" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img34.gif" - ALT="$ L_{max}$"> is the total number + functions contributing most to the reduction of \f$ S\f$ are chosen first + (TMultiDimFit::TestFunction). + + In case \f$ \phi \f$ isn't defined, an alternative method of + performing this second test is used: The \f$L^{\mbox{th}}\f$ + function \f$\textbf{f}_L\f$ is accepted if (refer also to equation (13)) + \f[ + \Delta S_L > \frac{S_{L-1}}{L_{max}-L} + \f] + where \f$ S_{L-1}\f$ is the sum of the \f$ L-1\f$ first residuals from the + \f$ L-1\f$ functions previously accepted; and \f$ L_{max}\f$ is the total number of functions allowed in the final expression of the fit (defined by user). - <P> - >From this we see, that by restricting <IMG - WIDTH="41" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img34.gif" - ALT="$ L_{max}$"> -- the number of + From this we see, that by restricting \f$ L_{max}\f$ -- the number of terms in the final model -- the fit is more difficult to perform, since the above selection criteria is more limiting. - <P> The more coefficients we evaluate, the more the sum of squares of - residuals <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> will be reduced. We can evaluate <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$"> before inverting - <!-- MATH - $\mathsf{B}$ - --> - <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img90.gif" - ALT="$ \mathsf{B}$"> as shown below. - - <P> - - <H2><A NAME="SECTION00034000000000000000"> - Coefficients and Coefficient Errors</A> - </H2> - - <P> - Having found a parameterization, that is the <IMG - WIDTH="19" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img24.gif" - ALT="$ F_l$">'s and <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img23.gif" - ALT="$ L$">, that - minimizes <IMG - WIDTH="15" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img26.gif" - ALT="$ S$">, we still need to determine the coefficients - <IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img25.gif" - ALT="$ c_l$">. However, it's a feature of how we choose the significant - functions, that the evaluation of the <IMG - WIDTH="16" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img25.gif" - ALT="$ c_l$">'s becomes trivial - [<A - HREF="TMultiFimFit.html#wind72">5</A>]. To derive <!-- MATH - $\mathbf{c}$ - --> - <IMG - WIDTH="12" HEIGHT="13" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img91.gif" - ALT="$ \mathbf{c}$">, we first note that - equation (<A HREF="TMultiFimFit.html#eq:wj">4</A>) can be written as - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:FF"></A><!-- MATH - \begin{equation} - \mathsf{F} = \mathsf{W}\mathsf{B} - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="60" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img92.gif" - ALT="$\displaystyle \mathsf{F} = \mathsf{W}\mathsf{B}$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (15)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> + residuals \f$ S\f$ will be reduced. We can evaluate \f$ S\f$ before inverting + \f$\mathsf{B}\f$ as shown below. + + ## Coefficients and Coefficient Errors + Having found a parameterization, that is the \f$ F_l\f$'s and \f$ L\f$, that + minimizes \f$ S\f$, we still need to determine the coefficients + \f$ c_l\f$. However, it's a feature of how we choose the significant + functions, that the evaluation of the \f$ c_l\f$'s becomes trivial [5]. To derive + \f$\mathbf{c}\f$, we first note that + equation (4) can be written as + \f[ + \mathsf{F} = \mathsf{W}\mathsf{B} + \f] where - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:bij"></A><!-- MATH - \begin{equation} - b_{ij} = \left\{\begin{array}{rcl} - \frac{\mathbf{f}_j \bullet \mathbf{w}_i}{\mathbf{w}_i^2} - & \mbox{if} & i < j\\ - 1 & \mbox{if} & i = j\\ - 0 & \mbox{if} & i > j - \end{array}\right. - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="187" HEIGHT="79" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img93.gif" - ALT="$\displaystyle b_{ij} = \left\{\begin{array}{rcl} \frac{\mathbf{f}_j \bullet \ma... - ...f} & i < j\ 1 & \mbox{if} & i = j\ 0 & \mbox{if} & i > j \end{array}\right.$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (16)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - Consequently, <!-- MATH - $\mathsf{B}$ - --> - <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img90.gif" - ALT="$ \mathsf{B}$"> is an upper triangle matrix, which can be + \f{eqnarray*}{ + b_{ij} = \frac{\mathbf{f}_j \bullet \mathbf{w}_i}{\mathbf{w}_i^2} + & \mbox{if} & i < j\\ + 1 & \mbox{if} & i = j\\ + 0 & \mbox{if} & i > j + \f} + Consequently, \f$\mathsf{B}\f$ is an upper triangle matrix, which can be readily inverted. So we now evaluate - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:FFF"></A><!-- MATH - \begin{equation} - \mathsf{F}\mathsf{B}^{-1} = \mathsf{W} - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="77" HEIGHT="35" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img94.gif" - ALT="$\displaystyle \mathsf{F}\mathsf{B}^{-1} = \mathsf{W}$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (17)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - The model <!-- MATH - $\mathsf{W}\mathbf{a}$ - --> - <IMG - WIDTH="28" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img63.gif" - ALT="$ \mathsf{W}\mathbf{a}$"> can therefore be written as - <!-- MATH - \begin{displaymath} - (\mathsf{F}\mathsf{B}^{-1})\mathbf{a} = - \mathsf{F}(\mathsf{B}^{-1}\mathbf{a})\,. - \end{displaymath} - --> - <P></P><DIV ALIGN="CENTER"> - <IMG - WIDTH="148" HEIGHT="35" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img95.gif" - ALT="$\displaystyle (\mathsf{F}\mathsf{B}^{-1})\mathbf{a} = - \mathsf{F}(\mathsf{B}^{-1}\mathbf{a}) . - $"> - </DIV><P></P> - The original model <!-- MATH - $\mathsf{F}\mathbf{c}$ - --> - <IMG - WIDTH="21" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img96.gif" - ALT="$ \mathsf{F}\mathbf{c}$"> is therefore identical with + \f[ + \mathsf{F}\mathsf{B}^{-1} = \mathsf{W} + \f] + The model \f$\mathsf{W}\mathbf{a}\f$ can therefore be written as + \f$(\mathsf{F}\mathsf{B}^{-1})\mathbf{a} = \mathsf{F}(\mathsf{B}^{-1}\mathbf{a})\,.\f$ + + The original model \f$\mathsf{F}\mathbf{c}\f$ is therefore identical with this if - <P></P> - <DIV ALIGN="CENTER"><A NAME="eq:id:cond"></A><!-- MATH - \begin{equation} - \mathbf{c} = \left(\mathsf{B}^{-1}\mathbf{a}\right) = - \left[\mathbf{a}^T\left(\mathsf{B}^{-1}\right)^T\right]^T\,. - \end{equation} - --> - <TABLE CELLPADDING="0" WIDTH="100%" ALIGN="CENTER"> - <TR VALIGN="MIDDLE"> - <TD NOWRAP ALIGN="CENTER"><IMG - WIDTH="214" HEIGHT="51" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img97.gif" - ALT="$\displaystyle \mathbf{c} = \left(\mathsf{B}^{-1}\mathbf{a}\right) = \left[\mathbf{a}^T\left(\mathsf{B}^{-1}\right)^T\right]^T .$"></TD> - <TD NOWRAP WIDTH="10" ALIGN="RIGHT"> - (18)</TD></TR> - </TABLE></DIV> - <BR CLEAR="ALL"><P></P> - The reason we use <!-- MATH - $\left(\mathsf{B}^{-1}\right)^T$ - --> - <IMG - WIDTH="56" HEIGHT="42" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img98.gif" - ALT="$ \left(\mathsf{B}^{-1}\right)^T$"> rather then - <!-- MATH - $\mathsf{B}^{-1}$ - --> - <IMG - WIDTH="32" HEIGHT="16" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img99.gif" - ALT="$ \mathsf{B}^{-1}$"> is to save storage, since - <!-- MATH - $\left(\mathsf{B}^{-1}\right)^T$ - --> - <IMG - WIDTH="56" HEIGHT="42" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img98.gif" - ALT="$ \left(\mathsf{B}^{-1}\right)^T$"> can be stored in the same matrix as - <!-- MATH - $\mathsf{B}$ - --> - <IMG - WIDTH="15" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img90.gif" - ALT="$ \mathsf{B}$"> - (<A NAME="tex2html17" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:MakeCoefficients"><TT>TMultiDimFit::MakeCoefficients</TT></A>). The errors in - the coefficients is calculated by inverting the curvature matrix - of the non-orthogonal functions <IMG - WIDTH="23" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img100.gif" - ALT="$ f_{lj}$"> [<A - HREF="TMultiFimFit.html#bevington">1</A>] - (<A NAME="tex2html18" - HREF=" - - ./TMultiDimFit.html#TMultiDimFit:MakeCoefficientErrors"><TT>TMultiDimFit::MakeCoefficientErrors</TT></A>). - - <P> - - <H2><A NAME="SECTION00035000000000000000"></A> - <A NAME="sec:considerations"></A><BR> - Considerations - </H2> - - <P> + \f[ + \mathbf{c} = \left(\mathsf{B}^{-1}\mathbf{a}\right) = + \left[\mathbf{a}^T\left(\mathsf{B}^{-1}\right)^T\right]^T\,. + \f] + The reason we use \f$\left(\mathsf{B}^{-1}\right)^T\f$ rather then + \f$\mathsf{B}^{-1}\f$ is to save storage, since \f$\left(\mathsf{B}^{-1}\right)^T\f$ + can be stored in the same matrix as \f$\mathsf{B}\f$ (TMultiDimFit::MakeCoefficients). + The errors in the coefficients is calculated by inverting the curvature matrix + of the non-orthogonal functions \f$ f_{lj}\f$ [1] (TMultiDimFit::MakeCoefficientErrors). + + ## Considerations It's important to realize that the training sample should be representive of the problem at hand, in particular along the borders of the region of interest. This is because the algorithm presented - here, is a <I>interpolation</I>, rahter then a <I>extrapolation</I> - [<A - HREF="TMultiFimFit.html#wind72">5</A>]. - - <P> - Also, the independent variables <IMG - WIDTH="18" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img101.gif" - ALT="$ x_{i}$"> need to be linear + here, is a *interpolation*, rahter then a *extrapolation* [5]. + + Also, the independent variables \f$ x_{i}\f$ need to be linear independent, since the procedure will perform poorly if they are not - [<A - HREF="TMultiFimFit.html#wind72">5</A>]. One can find an linear transformation from ones - original variables <IMG - WIDTH="16" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img102.gif" - ALT="$ \xi_{i}$"> to a set of linear independent variables - <IMG - WIDTH="18" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img101.gif" - ALT="$ x_{i}$">, using a <I>Principal Components Analysis</I> - <A NAME="tex2html19" - HREF="./TPrincipal.html">(see <TT>TPrincipal</TT>)</A>, and - then use the transformed variable as input to this class [<A - HREF="TMultiFimFit.html#wind72">5</A>] - [<A - HREF="TMultiFimFit.html#wind81">6</A>]. - - <P> + [5]. One can find an linear transformation from ones + original variables \f$ \xi_{i}\f$ to a set of linear independent variables + \f$ x_{i}\f$, using a *Principal Components Analysis* (see TPrincipal), and + then use the transformed variable as input to this class [5] [6]. + H. Wind also outlines a method for parameterising a multidimensional dependence over a multidimensional set of variables. An example - of the method from [<A - HREF="TMultiFimFit.html#wind72">5</A>], is a follows (please refer to - [<A - HREF="TMultiFimFit.html#wind72">5</A>] for a full discussion): - - <P> - - <OL> - <LI>Define <!-- MATH - $\mathbf{P} = (P_1, \ldots, P_5)$ - --> - <IMG - WIDTH="123" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img103.gif" - ALT="$ \mathbf{P} = (P_1, \ldots, P_5)$"> are the 5 dependent + of the method from [5], is a follows (please refer to + [5] for a full discussion): + + 1. Define \f$\mathbf{P} = (P_1, \ldots, P_5)\f$ are the 5 dependent quantities that define a track. - </LI> - <LI>Compute, for <IMG - WIDTH="21" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img10.gif" - ALT="$ M$"> different values of <!-- MATH - $\mathbf{P}$ - --> - <IMG - WIDTH="17" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img104.gif" - ALT="$ \mathbf{P}$">, the tracks + 2. Compute, for \f$ M\f$ different values of \f$\mathbf{P}\f$, the tracks through the magnetic field, and determine the corresponding - <!-- MATH - $\mathbf{x} = (x_1, \ldots, x_N)$ - --> - <IMG - WIDTH="123" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img105.gif" - ALT="$ \mathbf{x} = (x_1, \ldots, x_N)$">. - </LI> - <LI>Use the simulated observations to determine, with a simple - approximation, the values of <!-- MATH - $\mathbf{P}_j$ - --> - <IMG - WIDTH="23" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img106.gif" - ALT="$ \mathbf{P}_j$">. We call these values - <!-- MATH - $\mathbf{P}^\prime_j, j = 1, \ldots, M$ - --> - <IMG - WIDTH="122" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img107.gif" - ALT="$ \mathbf{P}^\prime_j, j = 1, \ldots, M$">. - </LI> - <LI>Determine from <!-- MATH - $\mathbf{x}$ - --> - <IMG - WIDTH="14" HEIGHT="13" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img9.gif" - ALT="$ \mathbf{x}$"> a set of at least five relevant - coordinates <!-- MATH - $\mathbf{x}^\prime$ - --> - <IMG - WIDTH="18" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img108.gif" - ALT="$ \mathbf{x}^\prime$">, using contrains, <I>or - alternative:</I> - </LI> - <LI>Perform a Principal Component Analysis (using - <A NAME="tex2html20" - HREF="./TPrincipal.html"><TT>TPrincipal</TT></A>), and use - - to get a linear transformation - <!-- MATH - $\mathbf{x} \rightarrow \mathbf{x}^\prime$ - --> - <IMG - WIDTH="53" HEIGHT="16" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img109.gif" - ALT="$ \mathbf{x} \rightarrow \mathbf{x}^\prime$">, so that - <!-- MATH - $\mathbf{x}^\prime$ - --> - <IMG - WIDTH="18" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img108.gif" - ALT="$ \mathbf{x}^\prime$"> are constrained and linear independent. - </LI> - <LI>Perform a Principal Component Analysis on - <!-- MATH - $Q_i = P_i / P^\prime_i\, i = 1, \ldots, 5$ - --> - <IMG - WIDTH="210" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img110.gif" - ALT="$ Q_i = P_i / P^prime_i i = 1, \ldots, 5$">, to get linear - indenpendent (among themselves, but not independent of - <!-- MATH - $\mathbf{x}$ - --> - <IMG - WIDTH="14" HEIGHT="13" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img9.gif" - ALT="$ \mathbf{x}$">) quantities <!-- MATH - $\mathbf{Q}^\prime$ - --> - <IMG - WIDTH="22" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img111.gif" - ALT="$ \mathbf{Q}^\prime$"> - </LI> - <LI>For each component <!-- MATH - $Q^\prime_i$ - --> - <IMG - WIDTH="22" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img112.gif" - ALT="$ Q^\prime_i$"> make a mutlidimensional fit, - using <!-- MATH - $\mathbf{x}^\prime$ - --> - <IMG - WIDTH="18" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img108.gif" - ALT="$ \mathbf{x}^\prime$"> as the variables, thus determing a set of - coefficents <!-- MATH - $\mathbf{c}_i$ - --> - <IMG - WIDTH="17" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img113.gif" - ALT="$ \mathbf{c}_i$">. - </LI> - </OL> - - <P> - To process data, using this parameterisation, do + \f$\mathbf{x} = (x_1, \ldots, x_N)\f$. + 3. Use the simulated observations to determine, with a simple + approximation, the values of \f$\mathbf{P}_j\f$. We call these values + \f$\mathbf{P}^\prime_j, j = 1, \ldots, M\f$. + 4. Determine from \f$\mathbf{x}\f$ a set of at least five relevant + coordinates \f$\mathbf{x}^\prime\f$, using contrains, *or + alternative:* + 5. Perform a Principal Component Analysis (using TPrincipal), and use + to get a linear transformation \f$\mathbf{x} \rightarrow \mathbf{x}^\prime\f$, so that + \f$\mathbf{x}^\prime\f$ are constrained and linear independent. + 6. Perform a Principal Component Analysis on + \f$Q_i = P_i / P^\prime_i\, i = 1, \ldots, 5\f$, to get linear + indenpendent (among themselves, but not independent of \f$\mathbf{x}\f$) quantities + \f$\mathbf{Q}^\prime\f$ + 7. For each component \f$Q^\prime_i\f$ make a mutlidimensional fit, + using \f$\mathbf{x}^\prime\f$ as the variables, thus determing a set of + coefficents \f$\mathbf{c}_i\f$. - <OL> - <LI>Test wether the observation <!-- MATH - $\mathbf{x}$ - --> - <IMG - WIDTH="14" HEIGHT="13" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img9.gif" - ALT="$ \mathbf{x}$"> within the domain of + To process data, using this parameterisation, do + 1. Test wether the observation \f$\mathbf{x}\f$ within the domain of the parameterization, using the result from the Principal Component Analysis. - </LI> - <LI>Determine <!-- MATH - $\mathbf{P}^\prime$ - --> - <IMG - WIDTH="21" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img114.gif" - ALT="$ \mathbf{P}^\prime$"> as before. - </LI> - <LI>Detetmine <!-- MATH - $\mathbf{x}^\prime$ - --> - <IMG - WIDTH="18" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img108.gif" - ALT="$ \mathbf{x}^\prime$"> as before. - </LI> - <LI>Use the result of the fit to determind <!-- MATH - $\mathbf{Q}^\prime$ - --> - <IMG - WIDTH="22" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img111.gif" - ALT="$ \mathbf{Q}^\prime$">. - </LI> - <LI>Transform back to <!-- MATH - $\mathbf{P}$ - --> - <IMG - WIDTH="17" HEIGHT="14" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img104.gif" - ALT="$ \mathbf{P}$"> from <!-- MATH - $\mathbf{Q}^\prime$ - --> - <IMG - WIDTH="22" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img111.gif" - ALT="$ \mathbf{Q}^\prime$">, using + 2. Determine \f$\mathbf{P}^\prime\f$ as before. + 3. Detetmine \f$\mathbf{x}^\prime\f$ as before. + 4. Use the result of the fit to determind \f$\mathbf{Q}^\prime\f$. + 5. Transform back to \f$\mathbf{P}\f$ from \f$\mathbf{Q}^\prime\f$, using the result from the Principal Component Analysis. - </LI> - </OL> - <P> - - <H2><A NAME="SECTION00036000000000000000"></A> - <A NAME="sec:testing"></A><BR> - Testing the parameterization - </H2> - - <P> + ## Testing the parameterization The class also provides functionality for testing the, over the - training sample, found parameterization - (<A NAME="tex2html21" - HREF=" - ./TMultiDimFit.html#TMultiDimFit:Fit"><TT>TMultiDimFit::Fit</TT></A>). This is done by passing - the class a test sample of <IMG - WIDTH="25" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img115.gif" - ALT="$ M_t$"> tuples of the form <!-- MATH - $(\mathbf{x}_{t,j}, - D_{t,j}, E_{t,j})$ - --> - <IMG - WIDTH="111" HEIGHT="31" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img116.gif" - ALT="$ (\mathbf{x}_{t,j}, - D_{t,j}, E_{t,j})$">, where <!-- MATH - $\mathbf{x}_{t,j}$ - --> - <IMG - WIDTH="29" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img117.gif" - ALT="$ \mathbf{x}_{t,j}$"> are the independent - variables, <IMG - WIDTH="33" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img118.gif" - ALT="$ D_{t,j}$"> the known, dependent quantity, and <IMG - WIDTH="31" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img119.gif" - ALT="$ E_{t,j}$"> is - the square error in <IMG - WIDTH="33" HEIGHT="29" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img118.gif" - ALT="$ D_{t,j}$"> - (<A NAME="tex2html22" - HREF=" - - ./TMultiDimFit.html#TMultiDimFit:AddTestRow"><TT>TMultiDimFit::AddTestRow</TT></A>). - - <P> - The parameterization is then evaluated at every <!-- MATH - $\mathbf{x}_t$ - --> - <IMG - WIDTH="19" HEIGHT="28" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img120.gif" - ALT="$ \mathbf{x}_t$"> in the + training sample, found parameterization (TMultiDimFit::Fit). This is done by passing + the class a test sample of \f$ M_\f$ tuples of the form + \f$(\mathbf{x}_{t,j},D_{t,j}, E_{t,j})\f$, where \f$\mathbf{x}_{t,j}\f$ are the independent + variables, \f$ D_{t,j}\f$ the known, dependent quantity, and \f$ E_{t,j}\f$ is + the square error in \f$ D_{t,j}\f$ (TMultiDimFit::AddTestRow). + + The parameterization is then evaluated at every \f$\mathbf{x}_t\f$ in the test sample, and - <!-- MATH - \begin{displaymath} - S_t \equiv \sum_{j=1}^{M_t} \left(D_{t,j} - - D_p\left(\mathbf{x}_{t,j}\right)\right)^2 - \end{displaymath} - --> - <P></P><DIV ALIGN="CENTER"> - <IMG - WIDTH="194" HEIGHT="66" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img121.gif" - ALT="$\displaystyle S_t \equiv \sum_{j=1}^{M_t} \left(D_{t,j} - - D_p\left(\mathbf{x}_{t,j}\right)\right)^2 - $"> - </DIV><P></P> + \f[ + S_t \equiv \sum_{j=1}^{M_t} \left(D_{t,j} - + D_p\left(\mathbf{x}_{t,j}\right)\right)^2 + \f] is evaluated. The relative error over the test sample - <!-- MATH - \begin{displaymath} - R_t = \frac{S_t}{\sum_{j=1}^{M_t} D_{t,j}^2} - \end{displaymath} - --> - <P></P><DIV ALIGN="CENTER"> - <IMG - WIDTH="118" HEIGHT="51" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img122.gif" - ALT="$\displaystyle R_t = \frac{S_t}{\sum_{j=1}^{M_t} D_{t,j}^2} - $"> - </DIV><P></P> - should not be to low or high compared to <IMG - WIDTH="16" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/multidimfit_img123.gif" - ALT="$ R$"> from the training + \f[ + R_t = \frac{S_t}{\sum_{j=1}^{M_t} D_{t,j}^2} + \f] + should not be to low or high compared to \f$ R\f$ from the training sample. Also, multiple correlation coefficient from both samples should be fairly close, otherwise one of the samples is not representive of - the problem. A large difference in the reduced <IMG - WIDTH="21" HEIGHT="33" ALIGN="MIDDLE" BORDER="0" - SRC="gif/multidimfit_img124.gif" - ALT="$ \chi^2$"> over the two + the problem. A large difference in the reduced \f$ \chi^2\f$ over the two samples indicate an over fit, and the maximum number of terms in the parameterisation should be reduced. - <P> - It's possible to use <A NAME="tex2html23" - HREF="./TMinuit.html"><I>Minuit</I></A> - [<A - HREF="TMultiFimFit.html#minuit">4</A>] to further improve the fit, using the test sample. + It's possible to use [4] to further improve the fit, using the test sample. - <P> - <DIV ALIGN="RIGHT"> Christian Holm - <BR> November 2000, NBI - - </DIV> - - <P> - - <H2><A NAME="SECTION00040000000000000000"> - Bibliography</A> - </H2><DL COMPACT><DD><P></P><DT><A NAME="bevington">1</A> - <DD> - Philip R. Bevington and D. Keith Robinson. - <BR><EM>Data Reduction and Error Analysis for the Physical Sciences</EM>. - <BR>McGraw-Hill, 2 edition, 1992. - - <P></P><DT><A NAME="mudifi">2</A> - <DD> - René Brun et al. - <BR>Mudifi. - <BR>Long writeup DD/75-23, CERN, 1980. - - <P></P><DT><A NAME="golub">3</A> - <DD> - Gene H. Golub and Charles F. van Loan. - <BR><EM>Matrix Computations</EM>. - <BR>John Hopkins Univeristy Press, Baltimore, 3 edition, 1996. - - <P></P><DT><A NAME="minuit">4</A> - <DD> - F. James. - <BR>Minuit. - <BR>Long writeup D506, CERN, 1998. - - <P></P><DT><A NAME="wind72">5</A> - <DD> - H. Wind. - <BR>Function parameterization. - <BR>In <EM>Proceedings of the 1972 CERN Computing and Data Processing - School</EM>, volume 72-21 of <EM>Yellow report</EM>. CERN, 1972. - - <P></P><DT><A NAME="wind81">6</A> - <DD> - H. Wind. - <BR>1. principal component analysis, 2. pattern recognition for track - finding, 3. interpolation and functional representation. - <BR>Yellow report EP/81-12, CERN, 1981. - </DL> - <pre> - */ -//End_Html -// + + ## Bibliography + - Philip R. Bevington and D. Keith Robinson. *Data Reduction and Error Analysis for + the Physical Sciences*. McGraw-Hill, 2 edition, 1992. + - R. Brun et al. *Long writeup DD/75-23*, CERN, 1980. + - Gene H. Golub and Charles F. van Loan. *Matrix Computations*. + John Hopkins Univeristy Press, Baltimore, 3 edition, 1996. + - F. James. *Minuit*. Long writeup D506, CERN, 1998. + - H. Wind. *Function parameterization*. Proceedings of the 1972 CERN Computing and Data Processing + School, volume 72-21 of Yellow report. CERN, 1972. + - H. Wind. 1. principal component analysis, 2. pattern recognition for track + finding, 3. interpolation and functional representation. Yellow report EP/81-12, CERN, 1981. + +[1]: TMultiFimFit.html#bevington +[2]: TMultiFimFit.html#mudifi +[4]: TMultiFimFit.html#minuit +[5]: TMultiFimFit.html#wind72 +[6]: TMultiFimFit.html#wind81 +[9]: TMultiFimFit.html#eq:dS2 +*/ + #include "Riostream.h" #include "TMultiDimFit.h" diff --git a/hist/hist/src/TMultiGraph.cxx b/hist/hist/src/TMultiGraph.cxx index 43a622e7cca..b6eb5a8fcb1 100644 --- a/hist/hist/src/TMultiGraph.cxx +++ b/hist/hist/src/TMultiGraph.cxx @@ -39,9 +39,8 @@ ClassImp(TMultiGraph) //////////////////////////////////////////////////////////////////////////////// -/* Begin_Html -<center><h2>TMultiGraph class</h2></center> - +/** \class TMultiGraph + \ingroup Hist A TMultiGraph is a collection of TGraph (or derived) objects. It allows to manipulate a set of graphs as a single entity. In particular, when drawn, the X and Y axis ranges are automatically computed such as all the graphs @@ -52,10 +51,9 @@ will be visible. The TMultiGraph owns the objects in the list. <p> The drawing options are the same as for TGraph. -Like for TGraph, the painting is performed thanks to the -<a href="http://root.cern.ch/root/html/TGraphPainter.html">TGraphPainter</a> -class. All details about the various painting options are given in -<a href="http://root.cern.ch/root/html/TGraphPainter.html">this class</a>. +Like for TGraph, the painting is performed thanks to the TGraphPainter +class. All details about the various painting options are given in this class. + Example: <pre> TGraph *gr1 = new TGraph(... @@ -67,7 +65,7 @@ Example: </pre> A special option <tt>3D</tt> allows to draw the graphs in a 3D space. See the following example: -End_Html + Begin_Macro(source) { c0 = new TCanvas("c1","multigraph L3",200,10,700,500); @@ -100,7 +98,7 @@ Begin_Macro(source) return c0; } End_Macro -Begin_Html + <p> The number of graphs in a multigraph can be retrieve with: <pre> @@ -115,7 +113,7 @@ otherwise the graph will be drawn with the option specified in <tt>TMultiGraph::Draw</tt>. <p> The following example shows how to fit a TMultiGraph. -End_Html + Begin_Macro(source) { TCanvas *c1 = new TCanvas("c1","c1",600,400); @@ -150,7 +148,7 @@ Begin_Macro(source) return c1; } End_Macro -Begin_Html + <p> The axis titles can be modified the following way: <p> @@ -166,18 +164,16 @@ The axis titles can be modified the following way: When the graphs in a TMultiGraph are fitted, the fit parameters boxes overlap. The following example shows how to make them all visible. -End_Html + Begin_Macro(source) ../../../tutorials/graphs/multigraph.C End_Macro -Begin_Html <p> The axis limits can be changed the like for TGraph. The same methods apply on the multigraph. Note the two differents ways to change limits on X and Y axis. -End_Html Begin_Macro(source) { TCanvas *c2 = new TCanvas("c2","c2",600,400); @@ -206,13 +202,11 @@ Begin_Macro(source) return c2; } End_Macro -Begin_Html + <p> -The method <a href="http://root.cern.ch/root/html/TPad.html#TPad:BuildLegend"> -<tt>TPad::BuildLegend</tt></a> is able to extract the graphs inside a +The method TPad::BuildLegend is able to extract the graphs inside a multigraph. The following example demonstrate this. -End_Html Begin_Macro(source) { TCanvas *c3 = new TCanvas("c3","c3",600, 400); @@ -270,10 +264,7 @@ Begin_Macro(source) return c3; } End_Macro -Begin_Html - - -End_Html */ +*/ //////////////////////////////////////////////////////////////////////////////// diff --git a/hist/hist/src/TPolyMarker.cxx b/hist/hist/src/TPolyMarker.cxx index e05b3deab29..887e95af7db 100644 --- a/hist/hist/src/TPolyMarker.cxx +++ b/hist/hist/src/TPolyMarker.cxx @@ -19,13 +19,13 @@ ClassImp(TPolyMarker) -//______________________________________________________________________________ -// -// a PolyMarker is defined by an array on N points in a 2-D space. -// At each point x[i], y[i] a marker is drawn. -// Marker attributes are managed by TAttMarker. -// See TMarker for the list of possible marker types. -// +/** \class TPolyMarker + \ingroup Hist +A PolyMarker is defined by an array on N points in a 2-D space. +At each point x[i], y[i] a marker is drawn. +Marker attributes are managed by TAttMarker. +See TMarker for the list of possible marker types. +*/ //////////////////////////////////////////////////////////////////////////////// diff --git a/hist/hist/src/TPrincipal.cxx b/hist/hist/src/TPrincipal.cxx index 9b24f54734c..76bb5d42a01 100644 --- a/hist/hist/src/TPrincipal.cxx +++ b/hist/hist/src/TPrincipal.cxx @@ -9,30 +9,16 @@ * For the list of contributors see $ROOTSYS/README/CREDITS. * *************************************************************************/ -//____________________________________________________________________ -//Begin_Html <!-- -/* --> - </pre> -<H1><A NAME="SECTION00010000000000000000"></A> -<A NAME="sec:lintra"></A> -<BR> +/** \class TPrincipal + \ingroup Hist Principal Components Analysis (PCA) -</H1> -<P> The current implementation is based on the LINTRA package from CERNLIB by R. Brun, H. Hansroul, and J. Kubler. The class has been implemented by Christian Holm Christensen in August 2000. -<P> - -<H2><A NAME="SECTION00011000000000000000"></A> -<A NAME="sec:intro1"></A> -<BR> -Introduction -</H2> +## Introduction -<P> In many applications of various fields of research, the treatment of large amounts of data requires powerful techniques capable of rapid data reduction and analysis. Usually, the quantities most @@ -42,7 +28,6 @@ then useful to have a way of selecting an optimal set of variables necessary for the recognition process and reducing the dimensionality of the problem, resulting in an easier classification procedure. -<P> This paper describes the implementation of one such method of feature selection, namely the principal components analysis. This multidimensional technique is well known in the field of pattern @@ -50,373 +35,119 @@ recognition and and its use in Particle Physics has been documented elsewhere (cf. H. Wind, <I>Function Parameterization</I>, CERN 72-21). -<P> - -<H2><A NAME="SECTION00012000000000000000"></A> -<A NAME="sec:overview"></A> -<BR> -Overview -</H2> - -<P> +## Overview Suppose we have prototypes which are trajectories of particles, passing through a spectrometer. If one measures the passage of the particle at say 8 fixed planes, the trajectory is described by an 8-component vector: -<BR><P></P> -<DIV ALIGN="CENTER"> - -<!-- MATH - \begin{displaymath} -\mathbf{x} = \left(x_0, x_1, \ldots, x_7\right) -\end{displaymath} - --> - - -<IMG - WIDTH="145" HEIGHT="31" BORDER="0" - SRC="gif/principal_img1.gif" - ALT="\begin{displaymath} -\mathbf{x} = \left(x_0, x_1, \ldots, x_7\right) -\end{displaymath}"> -</DIV> -<BR CLEAR="ALL"> -<P></P> +\f[ + \mathbf{x} = \left(x_0, x_1, \ldots, x_7\right) +\f] in 8-dimensional pattern space. -<P> One proceeds by generating a a representative tracks sample and -building up the covariance matrix <IMG - WIDTH="16" HEIGHT="16" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img2.gif" - ALT="$\mathsf{C}$">. Its eigenvectors and +building up the covariance matrix \f$\mathsf{C}\f$. Its eigenvectors and eigenvalues are computed by standard methods, and thus a new basis is obtained for the original 8-dimensional space the expansion of the prototypes, -<BR><P></P> -<DIV ALIGN="CENTER"> - -<!-- MATH - \begin{displaymath} -\mathbf{x}_m = \sum^7_{i=0} a_{m_i} \mathbf{e}_i -\quad -\mbox{where} -\quad -a_{m_i} = \mathbf{x}^T\bullet\mathbf{e}_i -\end{displaymath} - --> - - -<IMG - WIDTH="295" HEIGHT="58" BORDER="0" - SRC="gif/principal_img3.gif" - ALT="\begin{displaymath} -\mathbf{x}_m = \sum^7_{i=0} a_{m_i} \mathbf{e}_i -\quad -\mbox{where} -\quad -a_{m_i} = \mathbf{x}^T\bullet\mathbf{e}_i -\end{displaymath}"> -</DIV> -<BR CLEAR="ALL"> -<P></P> - -<P> -allows the study of the behavior of the coefficients <IMG - WIDTH="31" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img4.gif" - ALT="$a_{m_i}$"> for all +\f[ + \mathbf{x}_m = \sum^7_{i=0} a_{m_i} \mathbf{e}_i + \quad + \mbox{where} + \quad + a_{m_i} = \mathbf{x}^T\bullet\mathbf{e}_i +\f] +allows the study of the behavior of the coefficients \f$a_{m_i}\f$ for all the tracks of the sample. The eigenvectors which are insignificant for the trajectory description in the expansion will have their -corresponding coefficients <IMG - WIDTH="31" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img4.gif" - ALT="$a_{m_i}$"> close to zero for all the +corresponding coefficients \f$a_{m_i}\f$ close to zero for all the prototypes. -<P> On one hand, a reduction of the dimensionality is then obtained by omitting these least significant vectors in the subsequent analysis. -<P> On the other hand, in the analysis of real data, these least significant variables(?) can be used for the pattern recognition problem of extracting the valid combinations of coordinates describing a true trajectory from the set of all possible wrong combinations. -<P> The program described here performs this principal components analysis on a sample of data provided by the user. It computes the covariance matrix, its eigenvalues ands corresponding eigenvectors and exhibits -the behavior of the principal components (<IMG - WIDTH="31" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img4.gif" - ALT="$a_{m_i}$">), thus providing +the behavior of the principal components \f$a_{m_i}\f$, thus providing to the user all the means of understanding their data. -<P> - -<H2><A NAME="SECTION00013000000000000000"></A> -<A NAME="sec:method"></A> -<BR> -Principal Components Method -</H2> - -<P> -Let's consider a sample of <IMG - WIDTH="23" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img5.gif" - ALT="$M$"> prototypes each being characterized by -<IMG - WIDTH="18" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img6.gif" - ALT="$P$"> variables -<!-- MATH - $x_0, x_1, \ldots, x_{P-1}$ - --> -<IMG - WIDTH="107" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img7.gif" - ALT="$x_0, x_1, \ldots, x_{P-1}$">. Each prototype is a point, or a -column vector, in a <IMG - WIDTH="18" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img6.gif" - ALT="$P$">-dimensional <I>pattern space</I>. -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -\mathbf{x} = \left[\begin{array}{c} +## Principal Components Method +Let's consider a sample of \f$M\f$ prototypes each being characterized by +\f$P\f$ variables \f$x_0, x_1, \ldots, x_{P-1}\f$. Each prototype is a point, or a +column vector, in a \f$P\f$-dimensional *Pattern space*. +\f[ + \mathbf{x} = \left[\begin{array}{c} x_0\\x_1\\\vdots\\x_{P-1}\end{array}\right]\,, -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><IMG - WIDTH="102" HEIGHT="102" BORDER="0" - SRC="gif/principal_img8.gif" - ALT="\begin{displaymath} -\mathbf{x} = \left[\begin{array}{c} -x_0\\ x_1\\ \vdots\\ x_{P-1}\end{array}\right]\,, -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(1)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> -where each <IMG - WIDTH="23" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img9.gif" - ALT="$x_n$"> represents the particular value associated with the -<IMG - WIDTH="15" HEIGHT="16" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img10.gif" - ALT="$n$">-dimension. +\f] +where each \f$x_n\f$ represents the particular value associated with the +\f$n\f$-dimension. -<P> -Those <IMG - WIDTH="18" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img6.gif" - ALT="$P$"> variables are the quantities accessible to the +Those \f$P\f$ variables are the quantities accessible to the experimentalist, but are not necessarily the most significant for the classification purpose. -<P> -The <I>Principal Components Method</I> consists of applying a -<I>linear</I> transformation to the original variables. This +The *Principal Components Method* consists of applying a +*linear* transformation to the original variables. This transformation is described by an orthogonal matrix and is equivalent to a rotation of the original pattern space into a new set of coordinate vectors, which hopefully provide easier feature identification and dimensionality reduction. -<P> Let's define the covariance matrix: -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -\mathsf{C} = \left\langle\mathbf{y}\mathbf{y}^T\right\rangle +\f[ + \mathsf{C} = \left\langle\mathbf{y}\mathbf{y}^T\right\rangle \quad\mbox{where}\quad \mathbf{y} = \mathbf{x} - \left\langle\mathbf{x}\right\rangle\,, -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:C"></A><IMG - WIDTH="267" HEIGHT="37" BORDER="0" - SRC="gif/principal_img11.gif" - ALT="\begin{displaymath} -\mathsf{C} = \left\langle\mathbf{y}\mathbf{y}^T\right\rangl... -...athbf{y} = \mathbf{x} - \left\langle\mathbf{x}\right\rangle\,, -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(2)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> -and the brackets indicate mean value over the sample of <IMG - WIDTH="23" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img5.gif" - ALT="$M$"> +\f] +and the brackets indicate mean value over the sample of \f$M\f$ prototypes. -<P> -This matrix <IMG - WIDTH="16" HEIGHT="16" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img2.gif" - ALT="$\mathsf{C}$"> is real, positive definite, symmetric, and will +This matrix \f$\mathsf{C}\f$ is real, positive definite, symmetric, and will have all its eigenvalues greater then zero. It will now be show that among the family of all the complete orthonormal bases of the pattern space, the base formed by the eigenvectors of the covariance matrix and belonging to the largest eigenvalues, corresponds to the most significant features of the description of the original prototypes. -<P> -let the prototypes be expanded on into a set of <IMG - WIDTH="20" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img12.gif" - ALT="$N$"> basis vectors - -<!-- MATH - $\mathbf{e}_n, n=0,\ldots,N,N+1, \ldots, P-1$ - --> -<IMG - WIDTH="233" HEIGHT="32" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img13.gif" - ALT="$\mathbf{e}_n, n=0,\ldots,N,N+1, \ldots, P-1$">, -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -\mathbf{y}_i = \sum^N_{i=0} a_{i_n} \mathbf{e}_n, +let the prototypes be expanded on into a set of \f$N\f$ basis vectors +\f$\mathbf{e}_n, n=0,\ldots,N,N+1, \ldots, P-1\f$ +\f[ + \mathbf{y}_i = \sum^N_{i=0} a_{i_n} \mathbf{e}_n, \quad i = 1, \ldots, M, \quad N < P-1 -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:yi"></A><IMG - WIDTH="303" HEIGHT="58" BORDER="0" - SRC="gif/principal_img14.gif" - ALT="\begin{displaymath} -\mathbf{y}_i = \sum^N_{i=0} a_{i_n} \mathbf{e}_n, -\quad -i = 0, \ldots, M, -\quad -N < P-1 -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(3)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> - -<P> -The `best' feature coordinates <IMG - WIDTH="23" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img15.gif" - ALT="$\mathbf{e}_n$">, spanning a <I>feature - space</I>, will be obtained by minimizing the error due to this +\f] +The `best' feature coordinates \f$\mathbf{e}_n\f$, spanning a *feature +space*, will be obtained by minimizing the error due to this truncated expansion, i.e., -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -\min\left(E_N\right) = +\f[ + \min\left(E_N\right) = \min\left[\left\langle\left(\mathbf{y}_i - \sum^N_{i=0} a_{i_n} \mathbf{e}_n\right)^2\right\rangle\right] -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:mini"></A><IMG - WIDTH="306" HEIGHT="65" BORDER="0" - SRC="gif/principal_img16.gif" - ALT="\begin{displaymath} -\min\left(E_N\right) = -\min\left[\left\langle\left(\mathb... -...\sum^N_{i=0} a_{i_n} \mathbf{e}_n\right)^2\right\rangle\right] -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(4)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> +\f] with the conditions: -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -\mathbf{e}_k\bullet\mathbf{e}_j = \delta_{jk} = +\f[ + \mathbf{e}_k\bullet\mathbf{e}_j = \delta_{jk} = \left\{\begin{array}{rcl} 1 & \mbox{for} & k = j\\ 0 & \mbox{for} & k \neq j \end{array}\right. -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:ortocond"></A><IMG - WIDTH="240" HEIGHT="54" BORDER="0" - SRC="gif/principal_img17.gif" - ALT="\begin{displaymath} -\mathbf{e}_k\bullet\mathbf{e}_j = \delta_{jk} = -\left\{\b... -...for} & k = j\\ -0 & \mbox{for} & k \neq j -\end{array}\right. -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(5)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> - -<P> -Multiplying (<A HREF="prin_node1.html#eq:yi">3</A>) by -<!-- MATH - $\mathbf{e}^T_n$ - --> -<IMG - WIDTH="24" HEIGHT="38" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img18.gif" - ALT="$\mathbf{e}^T_n$"> using (<A HREF="prin_node1.html#eq:ortocond">5</A>), +\f] +Multiplying (3) by \f$\mathbf{e}^T_n\f$ using (5), we get -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -a_{i_n} = \mathbf{y}_i^T\bullet\mathbf{e}_n\,, -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:ai"></A><IMG - WIDTH="108" HEIGHT="31" BORDER="0" - SRC="gif/principal_img19.gif" - ALT="\begin{displaymath} -a_{i_n} = \mathbf{y}_i^T\bullet\mathbf{e}_n\,, -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(6)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> +\f[ + a_{i_n} = \mathbf{y}_i^T\bullet\mathbf{e}_n\,, +\f] so the error becomes -<BR> -<DIV ALIGN="CENTER"><A NAME="eq:error"></A> - -<!-- MATH - \begin{eqnarray} -E_N &=& +\f{eqnarray*}{ + E_N &=& \left\langle\left[\sum_{n=N+1}^{P-1} a_{i_n}\mathbf{e}_n\right]^2\right\rangle\nonumber\\ &=& \left\langle\left[\sum_{n=N+1}^{P-1} \mathbf{y}_i^T\bullet\mathbf{e}_n\mathbf{e}_n\right]^2\right\rangle\nonumber\\ @@ -424,202 +155,37 @@ E_N &=& \left\langle\sum_{n=N+1}^{P-1} \mathbf{e}_n^T\mathbf{y}_i\mathbf{y}_i^T\mathbf{e}_n\right\rangle\nonumber\\ &=& \sum_{n=N+1}^{P-1} \mathbf{e}_n^T\mathsf{C}\mathbf{e}_n -\end{eqnarray} - --> - -<TABLE ALIGN="CENTER" CELLPADDING="0" WIDTH="100%"> -<TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"><IMG - WIDTH="30" HEIGHT="32" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img20.gif" - ALT="$\displaystyle E_N$"></TD> -<TD ALIGN="CENTER" NOWRAP><IMG - WIDTH="18" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img21.gif" - ALT="$\textstyle =$"></TD> -<TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="151" HEIGHT="80" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img22.gif" - ALT="$\displaystyle \left\langle\left[\sum_{n=N+1}^{P-1} a_{i_n}\mathbf{e}_n\right]^2\right\rangle$"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> - </TD></TR> -<TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"> </TD> -<TD ALIGN="CENTER" NOWRAP><IMG - WIDTH="18" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img21.gif" - ALT="$\textstyle =$"></TD> -<TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="184" HEIGHT="80" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img23.gif" - ALT="$\displaystyle \left\langle\left[\sum_{n=N+1}^{P-1} \mathbf{y}_i^T\bullet\mathbf{e}_n\mathbf{e}_n\right]^2\right\rangle$"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> - </TD></TR> -<TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"> </TD> -<TD ALIGN="CENTER" NOWRAP><IMG - WIDTH="18" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img21.gif" - ALT="$\textstyle =$"></TD> -<TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="156" HEIGHT="69" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img24.gif" - ALT="$\displaystyle \left\langle\sum_{n=N+1}^{P-1} \mathbf{e}_n^T\mathbf{y}_i\mathbf{y}_i^T\mathbf{e}_n\right\rangle$"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> - </TD></TR> -<TR VALIGN="MIDDLE"><TD NOWRAP ALIGN="RIGHT"> </TD> -<TD ALIGN="CENTER" NOWRAP><IMG - WIDTH="18" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img21.gif" - ALT="$\textstyle =$"></TD> -<TD ALIGN="LEFT" NOWRAP><IMG - WIDTH="104" HEIGHT="69" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img25.gif" - ALT="$\displaystyle \sum_{n=N+1}^{P-1} \mathbf{e}_n^T\mathsf{C}\mathbf{e}_n$"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(7)</TD></TR> -</TABLE></DIV> -<BR CLEAR="ALL"><P></P> - -<P> -The minimization of the sum in (<A HREF="prin_node1.html#eq:error">7</A>) is obtained when each -term -<!-- MATH - $\mathbf{e}_n^\mathsf{C}\mathbf{e}_n$ - --> -<IMG - WIDTH="41" HEIGHT="38" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img26.gif" - ALT="$\mathbf{e}_n^\mathsf{C}\mathbf{e}_n$"> is minimum, since <IMG - WIDTH="16" HEIGHT="16" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img2.gif" - ALT="$\mathsf{C}$"> is +\f} +The minimization of the sum in (7) is obtained when each +term \f$\mathbf{e}_n^\mathsf{C}\mathbf{e}_n\f$ is minimum, since \f$\mathsf{C}\f$ is positive definite. By the method of Lagrange multipliers, and the -condition (<A HREF="prin_node1.html#eq:ortocond">5</A>), we get - -<P> - -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -E_N = \sum^{P-1}_{n=N+1} \left(\mathbf{e}_n^T\mathsf{C}\mathbf{e}_n - +condition (5), we get +\f[ + E_N = \sum^{P-1}_{n=N+1} \left(\mathbf{e}_n^T\mathsf{C}\mathbf{e}_n - l_n\mathbf{e}_n^T\bullet\mathbf{e}_n + l_n\right) -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:minerror"></A><IMG - WIDTH="291" HEIGHT="60" BORDER="0" - SRC="gif/principal_img27.gif" - ALT="\begin{displaymath} -E_N = \sum^{P-1}_{n=N+1} \left(\mathbf{e}_n^T\mathsf{C}\mathbf{e}_n - -l_n\mathbf{e}_n^T\bullet\mathbf{e}_n + l_n\right) -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(8)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> -The minimum condition -<!-- MATH - $\frac{dE_N}{d\mathbf{e}^T_n} = 0$ - --> -<IMG - WIDTH="68" HEIGHT="40" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img28.gif" - ALT="$\frac{dE_N}{d\mathbf{e}^T_n} = 0$"> leads to the +\f] +The minimum condition \f$\frac{dE_N}{d\mathbf{e}^T_n} = 0\f$ leads to the equation -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -\mathsf{C}\mathbf{e}_n = l_n\mathbf{e}_n\,, -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:Ce"></A><IMG - WIDTH="91" HEIGHT="30" BORDER="0" - SRC="gif/principal_img29.gif" - ALT="\begin{displaymath} -\mathsf{C}\mathbf{e}_n = l_n\mathbf{e}_n\,, -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(9)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> -which shows that <IMG - WIDTH="23" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img15.gif" - ALT="$\mathbf{e}_n$"> is an eigenvector of the covariance -matrix <IMG - WIDTH="16" HEIGHT="16" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img2.gif" - ALT="$\mathsf{C}$"> with eigenvalue <IMG - WIDTH="19" HEIGHT="32" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img30.gif" - ALT="$l_n$">. The estimated minimum error is +\f[ + \mathsf{C}\mathbf{e}_n = l_n\mathbf{e}_n\,, +\f] +which shows that \f$\mathbf{e}_n\f$ is an eigenvector of the covariance +matrix \f$\mathsf{C}\f$ with eigenvalue \f$l_n\f$. The estimated minimum error is then given by -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -E_N \sim \sum^{P-1}_{n=N+1} \mathbf{e}_n^T\bullet l_n\mathbf{e}_n +\f[ + E_N \sim \sum^{P-1}_{n=N+1} \mathbf{e}_n^T\bullet l_n\mathbf{e}_n = \sum^{P-1}_{n=N+1} l_n\,, -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:esterror"></A><IMG - WIDTH="264" HEIGHT="60" BORDER="0" - SRC="gif/principal_img31.gif" - ALT="\begin{displaymath} -E_N \sim \sum^{P-1}_{n=N+1} \mathbf{e}_n^T\bullet l_n\mathbf{e}_n -= \sum^{P-1}_{n=N+1} l_n\,, -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(10)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> -where -<!-- MATH - $l_n,\,n=N+1,\ldots,P$ - --> -<IMG - WIDTH="161" HEIGHT="32" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img32.gif" - ALT="$l_n,\,n=N+1,\ldots,P-1$"> are the eigenvalues associated with the -omitted eigenvectors in the expansion (<A HREF="prin_node1.html#eq:yi">3</A>). Thus, by choosing -the <IMG - WIDTH="20" HEIGHT="15" ALIGN="BOTTOM" BORDER="0" - SRC="gif/principal_img12.gif" - ALT="$N$"> largest eigenvalues, and their associated eigenvectors, the -error <IMG - WIDTH="30" HEIGHT="32" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img33.gif" - ALT="$E_N$"> is minimized. +\f] +where \f$l_n,\,n=N+1,\ldots,P\f$ \f$l_n,\,n=N+1,\ldots,P-1\f$ are the eigenvalues associated with the +omitted eigenvectors in the expansion (3). Thus, by choosing +the \f$N\f$ largest eigenvalues, and their associated eigenvectors, the +error \f$E_N\f$ is minimized. -<P> The transformation matrix to go from the pattern space to the feature -space consists of the ordered eigenvectors - -<!-- MATH - $\mathbf{e}_1,\ldots,\mathbf{e}_P$ - --> -<IMG - WIDTH="80" HEIGHT="30" ALIGN="MIDDLE" BORDER="0" - SRC="gif/principal_img34.gif" - ALT="$\mathbf{e}_0,\ldots,\mathbf{e}_{P-1}$"> for its columns -<BR> -<DIV ALIGN="RIGHT"> - - -<!-- MATH - \begin{equation} -\mathsf{T} = \left[ +space consists of the ordered eigenvectors \f$\mathbf{e}_1,\ldots,\mathbf{e}_P\f$ +\f$\mathbf{e}_0,\ldots,\mathbf{e}_{P-1}\f$ for its columns +\f[ + \mathsf{T} = \left[ \begin{array}{cccc} \mathbf{e}_0 & \mathbf{e}_1 & @@ -633,35 +199,13 @@ space consists of the ordered eigenvectors \vdots & \vdots & \ddots & \vdots \\ \mathbf{e}_{0_{P-1}} & \mathbf{e}_{1_{P-1}} & \cdots & \mathbf{e}_{{P-1}_{P-1}}\\ \end{array}\right] -\end{equation} - --> - -<TABLE WIDTH="100%" ALIGN="CENTER"> -<TR VALIGN="MIDDLE"><TD ALIGN="CENTER" NOWRAP><A NAME="eq:trans"></A><IMG - WIDTH="378" HEIGHT="102" BORDER="0" - SRC="gif/principal_img35.gif" - ALT="\begin{displaymath} -\mathsf{T} = \left[ -\begin{array}{cccc} -\mathbf{e}_0 & -\... -...bf{e}_{1_{P-1}} & \cdots & \mathbf{e}_{{P-1}_{P-1}}\\ -\end{array}\right] -\end{displaymath}"></TD> -<TD WIDTH=10 ALIGN="RIGHT"> -(11)</TD></TR> -</TABLE> -<BR CLEAR="ALL"></DIV><P></P> +\f] This is an orthogonal transformation, or rotation, of the pattern space and feature selection results in ignoring certain coordinates in the transformed space. - <p> - <DIV ALIGN="RIGHT"> - Christian Holm<br> - August 2000, CERN - </DIV> -<!--*/ -// -->End_Html + +Christian Holm August 2000, CERN +*/ // $Id$ // $Date: 2006/05/24 14:55:26 $ diff --git a/hist/hist/src/TProfile.cxx b/hist/hist/src/TProfile.cxx index a1ff9d9aef7..7648a4f2842 100644 --- a/hist/hist/src/TProfile.cxx +++ b/hist/hist/src/TProfile.cxx @@ -24,62 +24,59 @@ Bool_t TProfile::fgApproximate = kFALSE; ClassImp(TProfile) -//______________________________________________________________________________ -// -// Profile histograms are used to display the mean -// value of Y and its error for each bin in X. The displayed error is by default the -// standard error on the mean (i.e. the standard deviation divided by the sqrt(n) ) -// Profile histograms are in many cases an -// elegant replacement of two-dimensional histograms : the inter-relation of two -// measured quantities X and Y can always be visualized by a two-dimensional -// histogram or scatter-plot; its representation on the line-printer is not particularly -// satisfactory, except for sparse data. If Y is an unknown (but single-valued) -// approximate function of X, this function is displayed by a profile histogram with -// much better precision than by a scatter-plot. -// -// The following formulae show the cumulated contents (capital letters) and the values -// displayed by the printing or plotting routines (small letters) of the elements for bin J. -// -// 2 -// H(J) = sum Y E(J) = sum Y -// l(J) = sum l L(J) = sum l -// h(J) = H(J)/L(J) mean of Y, -// s(J) = sqrt(E(J)/L(J)- h(J)**2) standard deviation of Y (e.g. RMS) -// e(J) = s(J)/sqrt(L(J)) standard error on the mean -// -// The displayed bin content for bin J of a TProfile is always h(J). The corresponding bin error is by default -// e(J). In case the option "s" is used (in the constructor or by calling TProfile::BuildOptions) -// the displayed error is s(J) -// -// In the special case where s(J) is zero (eg, case of 1 entry only in one bin) -// the bin error e(J) is computed from the average of the s(J) for all bins if -// the static function TProfile::Approximate has been called. -// This simple/crude approximation was suggested in order to keep the bin -// during a fit operation. But note that this approximation is not the default behaviour. -// See also TProfile::BuildOptions for other error options and more detailed explanations -// -// Example of a profile histogram with its graphics output -//{ -// TCanvas *c1 = new TCanvas("c1","Profile histogram example",200,10,700,500); -// hprof = new TProfile("hprof","Profile of pz versus px",100,-4,4,0,20); -// Float_t px, py, pz; -// for ( Int_t i=0; i<25000; i++) { -// gRandom->Rannor(px,py); -// pz = px*px + py*py; -// hprof->Fill(px,pz,1); -// } -// hprof->Draw(); -//} -//Begin_Html -/* - <img src="gif/profile.gif"> - */ -//End_Html -// +/** \class TProfile + \ingroup Hist + Profile Historam. + Profile histograms are used to display the mean + value of Y and its error for each bin in X. The displayed error is by default the + standard error on the mean (i.e. the standard deviation divided by the sqrt(n) ) + Profile histograms are in many cases an + elegant replacement of two-dimensional histograms : the inter-relation of two + measured quantities X and Y can always be visualized by a two-dimensional + histogram or scatter-plot; its representation on the line-printer is not particularly + satisfactory, except for sparse data. If Y is an unknown (but single-valued) + approximate function of X, this function is displayed by a profile histogram with + much better precision than by a scatter-plot. + + The following formulae show the cumulated contents (capital letters) and the values + displayed by the printing or plotting routines (small letters) of the elements for bin J. + + 2 + H(J) = sum Y E(J) = sum Y + l(J) = sum l L(J) = sum l + h(J) = H(J)/L(J) mean of Y, + s(J) = sqrt(E(J)/L(J)- h(J)**2) standard deviation of Y (e.g. RMS) + e(J) = s(J)/sqrt(L(J)) standard error on the mean + + The displayed bin content for bin J of a TProfile is always h(J). The corresponding bin error is by default + e(J). In case the option "s" is used (in the constructor or by calling TProfile::BuildOptions) + the displayed error is s(J) + + In the special case where s(J) is zero (eg, case of 1 entry only in one bin) + the bin error e(J) is computed from the average of the s(J) for all bins if + the static function TProfile::Approximate has been called. + This simple/crude approximation was suggested in order to keep the bin + during a fit operation. But note that this approximation is not the default behaviour. + See also TProfile::BuildOptions for other error options and more detailed explanations + + Example of a profile histogram with its graphics output +~~~{.cpp} +{ + TCanvas *c1 = new TCanvas("c1","Profile histogram example",200,10,700,500); + hprof = new TProfile("hprof","Profile of pz versus px",100,-4,4,0,20); + Float_t px, py, pz; + for ( Int_t i=0; i<25000; i++) { + gRandom->Rannor(px,py); + pz = px*px + py*py; + hprof->Fill(px,pz,1); + } + hprof->Draw(); +} +~~~ +*/ //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Default constructor for Profile histograms*-*-*-*-*-*-*-*-* -///*-* ========================================== +/// Default constructor for Profile histograms TProfile::TProfile() : TH1D() { @@ -87,33 +84,31 @@ TProfile::TProfile() : TH1D() } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Default destructor for Profile histograms*-*-*-*-*-*-*-*-* -///*-* ========================================= +/// Default destructor for Profile histograms + TProfile::~TProfile() { } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Normal Constructor for Profile histograms*-*-*-*-*-*-*-*-*-* -///*-* ========================================== -/// -/// The first five parameters are similar to TH1D::TH1D. -/// All values of y are accepted at filling time. -/// To fill a profile histogram, one must use TProfile::Fill function. +/// Normal Constructor for Profile histograms +/// The first five parameters are similar to TH1D::TH1D. +/// All values of y are accepted at filling time. +/// To fill a profile histogram, one must use TProfile::Fill function. /// -/// Note that when filling the profile histogram the function Fill -/// checks if the variable y is betyween fYmin and fYmax. -/// If a minimum or maximum value is set for the Y scale before filling, -/// then all values below ymin or above ymax will be discarded. -/// Setting the minimum or maximum value for the Y scale before filling -/// has the same effect as calling the special TProfile constructor below -/// where ymin and ymax are specified. +/// Note that when filling the profile histogram the function Fill +/// checks if the variable y is betyween fYmin and fYmax. +/// If a minimum or maximum value is set for the Y scale before filling, +/// then all values below ymin or above ymax will be discarded. +/// Setting the minimum or maximum value for the Y scale before filling +/// has the same effect as calling the special TProfile constructor below +/// where ymin and ymax are specified. /// -/// H(J) is printed as the channel contents. The errors displayed are s(J) if CHOPT='S' -/// (spread option), or e(J) if CHOPT=' ' (error on mean). +/// H(J) is printed as the channel contents. The errors displayed are s(J) if CHOPT='S' +/// (spread option), or e(J) if CHOPT=' ' (error on mean). /// -/// See TProfile::BuildOptions for explanation of parameters +/// See TProfile::BuildOptions for explanation of parameters /// /// see also comments in the TH1 base class constructors @@ -124,11 +119,8 @@ TProfile::TProfile(const char *name,const char *title,Int_t nbins,Double_t xlow, } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Constructor for Profile histograms with variable bin size*-*-*-*-* -///*-* ========================================================= -/// -/// See TProfile::BuildOptions for more explanations on errors -/// +/// Constructor for Profile histograms with variable bin size +/// See TProfile::BuildOptions for more explanations on errors /// see also comments in the TH1 base class constructors TProfile::TProfile(const char *name,const char *title,Int_t nbins,const Float_t *xbins,Option_t *option) @@ -138,11 +130,8 @@ TProfile::TProfile(const char *name,const char *title,Int_t nbins,const Float_t } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Constructor for Profile histograms with variable bin size*-*-*-*-* -///*-* ========================================================= -/// -/// See TProfile::BuildOptions for more explanations on errors -/// +/// Constructor for Profile histograms with variable bin size +/// See TProfile::BuildOptions for more explanations on errors /// see also comments in the TH1 base class constructors TProfile::TProfile(const char *name,const char *title,Int_t nbins,const Double_t *xbins,Option_t *option) @@ -152,10 +141,8 @@ TProfile::TProfile(const char *name,const char *title,Int_t nbins,const Double_t } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Constructor for Profile histograms with variable bin size*-*-*-*-* -///*-* ========================================================= -/// -/// See TProfile::BuildOptions for more explanations on errors +/// Constructor for Profile histograms with variable bin size +/// See TProfile::BuildOptions for more explanations on errors /// /// see also comments in the TH1 base class constructors @@ -166,14 +153,13 @@ TProfile::TProfile(const char *name,const char *title,Int_t nbins,const Double_t } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Constructor for Profile histograms with range in y*-*-*-*-*-* -///*-* ================================================== -/// The first five parameters are similar to TH1D::TH1D. -/// Only the values of Y between ylow and yup will be considered at filling time. -/// ylow and yup will also be the maximum and minimum values -/// on the y scale when drawing the profile. +/// Constructor for Profile histograms with range in y +/// The first five parameters are similar to TH1D::TH1D. +/// Only the values of Y between ylow and yup will be considered at filling time. +/// ylow and yup will also be the maximum and minimum values +/// on the y scale when drawing the profile. /// -/// See TProfile::BuildOptions for more explanations on errors +/// See TProfile::BuildOptions for more explanations on errors /// /// see also comments in the TH1 base class constructors @@ -185,53 +171,46 @@ TProfile::TProfile(const char *name,const char *title,Int_t nbins,Double_t xlow, //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Set Profile histogram structure and options*-*-*-*-*-*-*-*-* -///*-* =========================================== -/// ymin: minimum value allowed for y -/// ymax: maximum value allowed for y +/// Set Profile histogram structure and options +/// \param[in] ymin minimum value allowed for y +/// \param[in] ymax maximum value allowed for y /// if (ymin = ymax = 0) there are no limits on the allowed y values (ymin = -inf, ymax = +inf) -/// -/// option: this is the option for the computation of the y error of the profile ( TProfile::GetBinError ) -/// possible values for the options are: -/// -/// -/// ' ' (Default) the bin errors are the standard error on the mean of Y = S(Y)/SQRT(N) -/// where S(Y) is the standard deviation (RMS) of the Y data in the bin -/// and N is the number of bin entries (from TProfile::GetBinEntries(ibin) ) -/// (i.e the errors are the standard error on the bin content of the profile) -/// -/// 's' Errors are the standard deviation of Y, S(Y) -/// -/// 'i' Errors are S(Y)/SQRT(N) (standard error on the mean as in the default) -/// The only difference is only when the standard deviation in Y is zero. -/// In this case the error a standard deviation = 1/SQRT(12) is assumed and the error is -/// 1./SQRT(12*N). -/// This approximation assumes that the Y values are integer (e.g. ADC counts) -/// and have an implicit uncertainty of y +/- 0.5. With the assumption that the probability that y -/// takes any value between y-0.5 and y+0.5 is uniform, its standard error is 1/SQRT(12) -/// -/// 'g' Errors are 1./SQRT(W) where W is the sum of the weights for the bin J -/// W is obtained as from TProfile::GetBinEntries(ibin) -/// This errors corresponds to the standard deviation of weighted mean where each -/// measurement Y is uncorrelated and has an error sigma, which is expressed in the -/// weight used to fill the Profile: w = 1/sigma^2 -/// The resulting error in TProfile is then 1./SQRT( Sum(1./sigma^2) ) +/// \param[in] option this is the option for the computation of the y error of the profile ( TProfile::GetBinError ) +/// possible values for the options are: +/// - ' ' (Default) the bin errors are the standard error on the mean of Y = S(Y)/SQRT(N) +/// where S(Y) is the standard deviation (RMS) of the Y data in the bin +/// and N is the number of bin entries (from TProfile::GetBinEntries(ibin) ) +/// (i.e the errors are the standard error on the bin content of the profile) +/// - 's' Errors are the standard deviation of Y, S(Y) +/// - 'i' Errors are S(Y)/SQRT(N) (standard error on the mean as in the default) +/// The only difference is only when the standard deviation in Y is zero. +/// In this case the error a standard deviation = 1/SQRT(12) is assumed and the error is +/// 1./SQRT(12*N). +/// This approximation assumes that the Y values are integer (e.g. ADC counts) +/// and have an implicit uncertainty of y +/- 0.5. With the assumption that the probability that y +/// takes any value between y-0.5 and y+0.5 is uniform, its standard error is 1/SQRT(12) +/// - 'g' Errors are 1./SQRT(W) where W is the sum of the weights for the bin J +/// W is obtained as from TProfile::GetBinEntries(ibin) +/// This errors corresponds to the standard deviation of weighted mean where each +/// measurement Y is uncorrelated and has an error sigma, which is expressed in the +/// weight used to fill the Profile: w = 1/sigma^2 +/// The resulting error in TProfile is then 1./SQRT( Sum(1./sigma^2) ) void TProfile::BuildOptions(Double_t ymin, Double_t ymax, Option_t *option) { // - // In the case of Profile filled weights and with TProfile::Sumw2() called, - // STD(Y) is the standard deviation of the weighted sample Y and N is in this case the - // number of effective entries (TProfile::GetBinEffectiveEntries(ibin) ) + // In the case of Profile filled weights and with TProfile::Sumw2() called, + // STD(Y) is the standard deviation of the weighted sample Y and N is in this case the + // number of effective entries (TProfile::GetBinEffectiveEntries(ibin) ) // - // If a bin has N data points all with the same value Y (especially - // possible when dealing with integers), the spread in Y for that bin - // is zero, and the uncertainty assigned is also zero, and the bin is - // ignored in making subsequent fits. - // To avoid this problem one can use an approximation for the standard deviation S(Y), - // by using the average of all the S(Y) of the other Profile bins. To use this approximation - // one must call before TProfile::Approximate - // This approximayion applies only for the default and the 's' options + // If a bin has N data points all with the same value Y (especially + // possible when dealing with integers), the spread in Y for that bin + // is zero, and the uncertainty assigned is also zero, and the bin is + // ignored in making subsequent fits. + // To avoid this problem one can use an approximation for the standard deviation S(Y), + // by using the average of all the S(Y) of the other Profile bins. To use this approximation + // one must call before TProfile::Approximate + // This approximayion applies only for the default and the 's' options // SetErrorOption(option); @@ -283,14 +262,11 @@ Bool_t TProfile::Add(const TH1 *h1, Double_t c1) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Replace contents of this profile by the addition of h1 and h2*-*-* -///*-* ============================================================= -/// -/// this = c1*h1 + c2*h2 -/// -/// c1 and c2 are considered as weights applied to the two summed profiles. -/// The operation acts therefore like merging the two profiles with a weight c1 and c2 +/// Replace contents of this profile by the addition of h1 and h2 +/// this = c1*h1 + c2*h2 /// +/// c1 and c2 are considered as weights applied to the two summed profiles. +/// The operation acts therefore like merging the two profiles with a weight c1 and c2 Bool_t TProfile::Add(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2) { @@ -311,14 +287,14 @@ Bool_t TProfile::Add(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2) //////////////////////////////////////////////////////////////////////////////// -/// static function +/// static function /// set the fgApproximate flag. When the flag is true, the function GetBinError /// will approximate the bin error with the average profile error on all bins /// in the following situation only -/// - the number of bins in the profile is less than 1002 -/// - the bin number of entries is small ( <5) -/// - the estimated bin error is extremely small compared to the bin content -/// (see TProfile::GetBinError) +/// - the number of bins in the profile is less than 1002 +/// - the bin number of entries is small ( <5) +/// - the estimated bin error is extremely small compared to the bin content +/// (see TProfile::GetBinError) void TProfile::Approximate(Bool_t approx) { @@ -327,9 +303,9 @@ void TProfile::Approximate(Bool_t approx) //////////////////////////////////////////////////////////////////////////////// /// Fill histogram with all entries in the buffer. -/// action = -1 histogram is reset and refilled from the buffer (called by THistPainter::Paint) -/// action = 0 histogram is filled from the buffer -/// action = 1 histogram is filled and buffer is deleted +/// - action = -1 histogram is reset and refilled from the buffer (called by THistPainter::Paint) +/// - action = 0 histogram is filled from the buffer +/// - action = 1 histogram is filled and buffer is deleted /// The buffer is automatically deleted when the number of entries /// in the buffer is greater than the number of entries in the histogram @@ -385,10 +361,10 @@ Int_t TProfile::BufferEmpty(Int_t action) //////////////////////////////////////////////////////////////////////////////// /// accumulate arguments in buffer. When buffer is full, empty the buffer -/// fBuffer[0] = number of entries in buffer -/// fBuffer[1] = w of first entry -/// fBuffer[2] = x of first entry -/// fBuffer[3] = y of first entry +/// - fBuffer[0] = number of entries in buffer +/// - fBuffer[1] = w of first entry +/// - fBuffer[2] = x of first entry +/// - fBuffer[3] = y of first entry Int_t TProfile::BufferFill(Double_t x, Double_t y, Double_t w) { @@ -415,8 +391,7 @@ Int_t TProfile::BufferFill(Double_t x, Double_t y, Double_t w) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*Copy a Profile histogram to a new profile histogram*-*-*-*-* -///*-* =================================================== +/// Copy a Profile histogram to a new profile histogram void TProfile::Copy(TObject &obj) const { @@ -455,10 +430,8 @@ Bool_t TProfile::Divide(TF1 *, Double_t ) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Divide this profile by h1*-*-*-*-*-*-*-*-*-*-*-*-* -///*-* ========================= -/// -/// this = this/h1 +/// Divide this profile by h1 +/// this = this/h1 /// This function accepts to divide a TProfile by a histogram /// /// The function return kFALSE if the divide operation failed @@ -480,16 +453,16 @@ Bool_t TProfile::Divide(const TH1 *h1) Int_t nbinsx = GetNbinsX(); - //*-*- Check profile compatibility + //- Check profile compatibility if (nbinsx != p1->GetNbinsX()) { Error("Divide","Attempt to divide profiles with different number of bins"); return kFALSE; } - //*-*- Reset statistics + //- Reset statistics fEntries = fTsumw = fTsumw2 = fTsumwx = fTsumwx2 = fTsumwy = fTsumwy2 = 0; - //*-*- Loop on bins (including underflows/overflows) + //- Loop on bins (including underflows/overflows) Int_t bin; Double_t *cu1=0, *er1=0, *en1=0; Double_t e0,e1,c12; @@ -537,10 +510,8 @@ Bool_t TProfile::Divide(const TH1 *h1) //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Replace contents of this profile by the division of h1 by h2*-*-* -///*-* ============================================================ -/// -/// this = c1*h1/(c2*h2) +/// Replace contents of this profile by the division of h1 by h2 +/// this = c1*h1/(c2*h2) /// /// The function return kFALSE if the divide operation failed @@ -569,7 +540,7 @@ Bool_t TProfile::Divide(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2, if (fBuffer) BufferEmpty(1); Int_t nbinsx = GetNbinsX(); - //*-*- Check histogram compatibility + //- Check histogram compatibility if (nbinsx != p1->GetNbinsX() || nbinsx != p2->GetNbinsX()) { Error("Divide","Attempt to divide profiles with different number of bins"); return kFALSE; @@ -586,10 +557,10 @@ Bool_t TProfile::Divide(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2, printf(" TH1D *p2 = h2->ProjectionX();\n"); printf(" p1->Divide(p2);\n"); - //*-*- Reset statistics + //- Reset statistics fEntries = fTsumw = fTsumw2 = fTsumwx = fTsumwx2 = 0; - //*-*- Loop on bins (including underflows/overflows) + //- Loop on bins (including underflows/overflows) Int_t bin; Double_t *cu1 = p1->GetW(); Double_t *cu2 = p2->GetW(); @@ -645,8 +616,7 @@ Bool_t TProfile::Divide(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2, } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Fill a Profile histogram (no weights)*-*-*-*-*-*-*-* -///*-* ===================================== +/// Fill a Profile histogram (no weights) Int_t TProfile::Fill(Double_t x, Double_t y) { @@ -706,8 +676,7 @@ Int_t TProfile::Fill(const char *namex, Double_t y) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Fill a Profile histogram with weights*-*-*-*-*-*-*-* -///*-* ===================================== +/// Fill a Profile histogram with weights Int_t TProfile::Fill(Double_t x, Double_t y, Double_t w) { @@ -773,8 +742,7 @@ Int_t TProfile::Fill(const char *namex, Double_t y, Double_t w) //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Fill a Profile histogram with weights*-*-*-*-*-*-*-* -///*-* ===================================== +/// Fill a Profile histogram with weights void TProfile::FillN(Int_t ntimes, const Double_t *x, const Double_t *y, const Double_t *w, Int_t stride) { @@ -822,8 +790,7 @@ void TProfile::FillN(Int_t ntimes, const Double_t *x, const Double_t *y, const D } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Return bin content of a Profile histogram*-*-*-*-*-*-*-*-*-* -///*-* ========================================= +/// Return bin content of a Profile histogram Double_t TProfile::GetBinContent(Int_t bin) const { @@ -836,8 +803,7 @@ Double_t TProfile::GetBinContent(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Return bin entries of a Profile histogram*-*-*-*-*-*-*-*-*-* -///*-* ========================================= +/// Return bin entries of a Profile histogram Double_t TProfile::GetBinEntries(Int_t bin) const { @@ -848,14 +814,14 @@ Double_t TProfile::GetBinEntries(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -/// Return bin effective entries for a weighted filled Profile histogram. -/// In case of an unweighted profile, it is equivalent to the number of entries per bin -/// The effective entries is defined as the square of the sum of the weights divided by the -/// sum of the weights square. -/// TProfile::Sumw2() must be called before filling the profile with weights. -/// Only by calling this method the sum of the square of the weights per bin is stored. +/// Return bin effective entries for a weighted filled Profile histogram. +/// In case of an unweighted profile, it is equivalent to the number of entries per bin +/// The effective entries is defined as the square of the sum of the weights divided by the +/// sum of the weights square. +/// TProfile::Sumw2() must be called before filling the profile with weights. +/// Only by calling this method the sum of the square of the weights per bin is stored. /// -///*-* ========================================= +/// ========================================= Double_t TProfile::GetBinEffectiveEntries(Int_t bin) const { @@ -863,30 +829,29 @@ Double_t TProfile::GetBinEffectiveEntries(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -/// *-*-*-*-*-*-*Return bin error of a Profile histogram*-*-*-*-*-*-*-*-*-* -/// *-* ======================================= +/// Return bin error of a Profile histogram /// /// Computing errors: A moving field -/// ================================= +/// /// The computation of errors for a TProfile has evolved with the versions /// of ROOT. The difficulty is in computing errors for bins with low statistics. /// - prior to version 3.00, we had no special treatment of low statistic bins. -/// As a result, these bins had huge errors. The reason is that the -/// expression eprim2 is very close to 0 (rounding problems) or 0. +/// As a result, these bins had huge errors. The reason is that the +/// expression eprim2 is very close to 0 (rounding problems) or 0. /// - in version 3.00 (18 Dec 2000), the algorithm is protected for values of -/// eprim2 very small and the bin errors set to the average bin errors, following -/// recommendations from a group of users. +/// eprim2 very small and the bin errors set to the average bin errors, following +/// recommendations from a group of users. /// - in version 3.01 (19 Apr 2001), it is realized that the algorithm above -/// should be applied only to low statistic bins. +/// should be applied only to low statistic bins. /// - in version 3.02 (26 Sep 2001), the same group of users recommend instead -/// to take two times the average error on all bins for these low -/// statistics bins giving a very small value for eprim2. +/// to take two times the average error on all bins for these low +/// statistics bins giving a very small value for eprim2. /// - in version 3.04 (Nov 2002), the algorithm is modified/protected for the case -/// when a TProfile is projected (ProjectionX). The previous algorithm -/// generated a N^2 problem when projecting a TProfile with a large number of -/// bins (eg 100000). +/// when a TProfile is projected (ProjectionX). The previous algorithm +/// generated a N^2 problem when projecting a TProfile with a large number of +/// bins (eg 100000). /// - in version 3.05/06, a new static function TProfile::Approximate -/// is introduced to enable or disable (default) the approximation. +/// is introduced to enable or disable (default) the approximation. /// /// Ideas for improvements of this algorithm are welcome. No suggestions /// received since our call for advice to roottalk in Jul 2002. @@ -898,8 +863,7 @@ Double_t TProfile::GetBinError(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Return option to compute profile errors*-*-*-*-*-*-*-*-* -///*-* ======================================= +/// Return option to compute profile errors Option_t *TProfile::GetErrorOption() const { @@ -988,14 +952,14 @@ void TProfile::LabelsInflate(Option_t *options) } //////////////////////////////////////////////////////////////////////////////// -/// Set option(s) to draw axis with labels -/// option = "a" sort by alphabetic order -/// = ">" sort by decreasing values -/// = "<" sort by increasing values -/// = "h" draw labels horizonthal -/// = "v" draw labels vertical -/// = "u" draw labels up (end of label right adjusted) -/// = "d" draw labels down (start of label left adjusted) +/// Set option(s) to draw axis with labels +/// option = "a" sort by alphabetic order +/// = ">" sort by decreasing values +/// = "<" sort by increasing values +/// = "h" draw labels horizonthal +/// = "v" draw labels vertical +/// = "u" draw labels up (end of label right adjusted) +/// = "d" draw labels down (start of label left adjusted) void TProfile::LabelsOption(Option_t *option, Option_t * /*ax */) { @@ -1155,7 +1119,7 @@ Bool_t TProfile::Multiply(TF1 *f1, Double_t c1) Int_t nbinsx = GetNbinsX(); - //*-*- Add statistics + //- Add statistics Double_t xx[1], cf1, ac1 = TMath::Abs(c1); Double_t s1[10]; Int_t i; @@ -1165,7 +1129,7 @@ Bool_t TProfile::Multiply(TF1 *f1, Double_t c1) SetMinimum(); SetMaximum(); - //*-*- Loop on bins (including underflows/overflows) + //- Loop on bins (including underflows/overflows) Int_t bin; for (bin=0;bin<=nbinsx+1;bin++) { xx[0] = fXaxis.GetBinCenter(bin); @@ -1183,11 +1147,9 @@ Bool_t TProfile::Multiply(TF1 *f1, Double_t c1) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Multiply this profile by h1*-*-*-*-*-*-*-*-*-*-*-*-* -///*-* ============================= -/// -/// this = this*h1 +/// Multiply this profile by h1 /// +/// this = this*h1 Bool_t TProfile::Multiply(const TH1 *) { @@ -1197,11 +1159,9 @@ Bool_t TProfile::Multiply(const TH1 *) //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Replace contents of this profile by multiplication of h1 by h2*-* -///*-* ================================================================ -/// -/// this = (c1*h1)*(c2*h2) +/// Replace contents of this profile by multiplication of h1 by h2 /// +/// this = (c1*h1)*(c2*h2) Bool_t TProfile::Multiply(const TH1 *, const TH1 *, Double_t, Double_t, Option_t *) { @@ -1210,26 +1170,25 @@ Bool_t TProfile::Multiply(const TH1 *, const TH1 *, Double_t, Double_t, Option_t } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Project this profile into a 1-D histogram along X*-*-*-*-*-*-* -///*-* ================================================= +/// Project this profile into a 1-D histogram along X /// -/// The projection is always of the type TH1D. +/// The projection is always of the type TH1D. /// -/// if option "E" is specified the errors of the projected histogram are computed and set -/// to be equal to the errors of the profile. -/// Option "E" is defined as the default one in the header file. -/// if option "" is specified the histogram errors are simply the sqrt of its content -/// if option "B" is specified, the content of bin of the returned histogram -/// will be equal to the GetBinEntries(bin) of the profile, -/// otherwise (default) it will be equal to GetBinContent(bin) -/// if option "C=E" the bin contents of the projection are set to the -/// bin errors of the profile -/// if option "W" is specified the bin content of the projected histogram is set to the -/// product of the bin content of the profile and the entries. -/// With this option the returned histogram will be equivalent to the one obtained by -/// filling directly a TH1D using the 2-nd value as a weight. -/// This makes sense only for profile filled with weights =1. If not, the error of the -/// projected histogram obtained with this option will not be correct. +/// - if option "E" is specified the errors of the projected histogram are computed and set +/// to be equal to the errors of the profile. +/// Option "E" is defined as the default one in the header file. +/// - if option "" is specified the histogram errors are simply the sqrt of its content +/// - if option "B" is specified, the content of bin of the returned histogram +/// will be equal to the GetBinEntries(bin) of the profile, +/// otherwise (default) it will be equal to GetBinContent(bin) +/// - if option "C=E" the bin contents of the projection are set to the +/// bin errors of the profile +/// - if option "W" is specified the bin content of the projected histogram is set to the +/// product of the bin content of the profile and the entries. +/// With this option the returned histogram will be equivalent to the one obtained by +/// filling directly a TH1D using the 2-nd value as a weight. +/// This makes sense only for profile filled with weights =1. If not, the error of the +/// projected histogram obtained with this option will not be correct. TH1D *TProfile::ProjectionX(const char *name, Option_t *option) const { @@ -1317,40 +1276,40 @@ void TProfile::PutStats(Double_t *stats) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Rebin this profile grouping ngroup bins together*-*-*-*-*-*-*-*-* -///*-* ================================================ -/// -case 1 xbins=0 -/// if newname is not blank a new temporary profile hnew is created. -/// else the current profile is modified (default) -/// The parameter ngroup indicates how many bins of this have to me merged -/// into one bin of hnew -/// If the original profile has errors stored (via Sumw2), the resulting -/// profile has new errors correctly calculated. +/// Rebin this profile grouping ngroup bins together +/// ## case 1 xbins=0 +/// if newname is not blank a new temporary profile hnew is created. +/// else the current profile is modified (default) +/// The parameter ngroup indicates how many bins of this have to me merged +/// into one bin of hnew +/// If the original profile has errors stored (via Sumw2), the resulting +/// profile has new errors correctly calculated. +/// +/// examples: if hp is an existing TProfile histogram with 100 bins +/// hp->Rebin(); //merges two bins in one in hp: previous contents of hp are lost +/// hp->Rebin(5); //merges five bins in one in hp +/// TProfile *hnew = hp->Rebin(5,"hnew"); // creates a new profile hnew +/// //merging 5 bins of hp in one bin /// -/// examples: if hp is an existing TProfile histogram with 100 bins -/// hp->Rebin(); //merges two bins in one in hp: previous contents of hp are lost -/// hp->Rebin(5); //merges five bins in one in hp -/// TProfile *hnew = hp->Rebin(5,"hnew"); // creates a new profile hnew -/// //merging 5 bins of hp in one bin +/// NOTE: If ngroup is not an exact divider of the number of bins, +/// the top limit of the rebinned profile is changed +/// to the upper edge of the bin=newbins*ngroup and the corresponding +/// bins are added to the overflow bin. +/// Statistics will be recomputed from the new bin contents. /// -/// NOTE: If ngroup is not an exact divider of the number of bins, -/// the top limit of the rebinned profile is changed -/// to the upper edge of the bin=newbins*ngroup and the corresponding -/// bins are added to the overflow bin. -/// Statistics will be recomputed from the new bin contents. +/// ## case 2 xbins!=0 +/// a new profile is created (you should specify newname). +/// The parameter ngroup is the number of variable size bins in the created profile +/// The array xbins must contain ngroup+1 elements that represent the low-edge +/// of the bins. +/// The data of the old bins are added to the new bin which contains the bin center +/// of the old bins. It is possible that information from the old binning are attached +/// to the under-/overflow bins of the new binning. /// -/// -case 2 xbins!=0 -/// a new profile is created (you should specify newname). -/// The parameter ngroup is the number of variable size bins in the created profile -/// The array xbins must contain ngroup+1 elements that represent the low-edge -/// of the bins. -/// The data of the old bins are added to the new bin which contains the bin center -/// of the old bins. It is possible that information from the old binning are attached -/// to the under-/overflow bins of the new binning. +/// examples: if hp is an existing TProfile with 100 bins /// -/// examples: if hp is an existing TProfile with 100 bins -/// Double_t xbins[25] = {...} array of low-edges (xbins[25] is the upper edge of last bin -/// hp->Rebin(24,"hpnew",xbins); //creates a new variable bin size profile hpnew +/// Double_t xbins[25] = {...} array of low-edges (xbins[25] is the upper edge of last bin +/// hp->Rebin(24,"hpnew",xbins); //creates a new variable bin size profile hpnew TH1 *TProfile::Rebin(Int_t ngroup, const char*newname, const Double_t *xbins) { @@ -1375,7 +1334,7 @@ TH1 *TProfile::Rebin(Int_t ngroup, const char*newname, const Double_t *xbins) } else { // in the case of xbins given (rebinning in variable bins) ngroup is the new number of bins. - // and number of grouped bins is not constant. + // and number of grouped bins is not constant. // when looping for setting the contents for the new histogram we // need to loop on all bins of original histogram. Set then ngroup=nbins newbins = ngroup; @@ -1536,8 +1495,7 @@ void TProfile::ExtendAxis(Double_t x, TAxis *axis) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Reset contents of a Profile histogram*-*-*-*-*-*-*-*-* -///*-* ===================================== +/// Reset contents of a Profile histogram void TProfile::Reset(Option_t *option) { @@ -1628,10 +1586,9 @@ void TProfile::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/) } //////////////////////////////////////////////////////////////////////////////// -/// *-*-*-*-*Multiply this profile by a constant c1*-*-*-*-*-*-*-*-* -/// *-* ====================================== +/// Multiply this profile by a constant c1 /// -/// this = c1*this +/// this = c1*this /// /// This function uses the services of TProfile::Add /// @@ -1642,8 +1599,7 @@ void TProfile::Scale(Double_t c1, Option_t * option) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*Set the number of entries in bin*-*-*-*-*-*-*-*-*-*-*-* -///*-* ================================ +/// Set the number of entries in bin void TProfile::SetBinEntries(Int_t bin, Double_t w) { @@ -1651,8 +1607,7 @@ void TProfile::SetBinEntries(Int_t bin, Double_t w) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*Redefine x axis parameters*-*-*-*-*-*-*-*-*-*-*-* -///*-* =========================== +/// Redefine x axis parameters void TProfile::SetBins(Int_t nx, Double_t xmin, Double_t xmax) { @@ -1662,8 +1617,7 @@ void TProfile::SetBins(Int_t nx, Double_t xmin, Double_t xmax) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*Redefine x axis parameters*-*-*-*-*-*-*-*-*-*-*-* -///*-* =========================== +/// Redefine x axis parameters void TProfile::SetBins(Int_t nx, const Double_t *xbins) { @@ -1703,30 +1657,26 @@ void TProfile::SetBuffer(Int_t buffersize, Option_t *) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Set option to compute profile errors*-*-*-*-*-*-*-*-* -///*-* ===================================== -/// -/// The computation of the bin errors is based on the parameter option: -/// option: -/// ' ' (Default) The bin errors are the standard error on the mean of the bin profiled values (Y), -/// i.e. the standard error of the bin contents. -/// Note that if TProfile::Approximate() is called, an approximation is used when -/// the spread in Y is 0 and the number of bin entries is > 0 -/// -/// 's' The bin errors are the standard deviations of the Y bin values -/// Note that if TProfile::Approximate() is called, an approximation is used when -/// the spread in Y is 0 and the number of bin entries is > 0 -/// -/// 'i' Errors are as in default case (standard errors of the bin contents) -/// The only difference is for the case when the spread in Y is zero. -/// In this case for N > 0 the error is 1./SQRT(12.*N) +/// Set option to compute profile errors /// -/// 'g' Errors are 1./SQRT(W) for W not equal to 0 and 0 for W = 0. -/// W is the sum in the bin of the weights of the profile. -/// This option is for combining measurements y +/- dy, -/// and the profile is filled with values y and weights w = 1/dy**2 +/// The computation of the bin errors is based on the parameter option: +/// option: +/// -' ' (Default) The bin errors are the standard error on the mean of the bin profiled values (Y), +/// i.e. the standard error of the bin contents. +/// Note that if TProfile::Approximate() is called, an approximation is used when +/// the spread in Y is 0 and the number of bin entries is > 0 +/// -'s' The bin errors are the standard deviations of the Y bin values +/// Note that if TProfile::Approximate() is called, an approximation is used when +/// the spread in Y is 0 and the number of bin entries is > 0 +/// -'i' Errors are as in default case (standard errors of the bin contents) +/// The only difference is for the case when the spread in Y is zero. +/// In this case for N > 0 the error is 1./SQRT(12.*N) +/// -'g' Errors are 1./SQRT(W) for W not equal to 0 and 0 for W = 0. +/// W is the sum in the bin of the weights of the profile. +/// This option is for combining measurements y +/- dy, +/// and the profile is filled with values y and weights w = 1/dy**2 /// -/// See TProfile::BuildOptions for a detailed explanation of all options +/// See TProfile::BuildOptions for a detailed explanation of all options void TProfile::SetErrorOption(Option_t *option) { @@ -1767,13 +1717,13 @@ void TProfile::Streamer(TBuffer &R__b) } } //////////////////////////////////////////////////////////////////////////////// -/// Create/delete structure to store sum of squares of weights per bin *-*-*-*-*-*-*-* -/// This is needed to compute the correct statistical quantities -/// of a profile filled with weights +/// Create/delete structure to store sum of squares of weights per bin --- +/// This is needed to compute the correct statistical quantities +/// of a profile filled with weights /// /// -/// This function is automatically called when the histogram is created -/// if the static function TH1::SetDefaultSumw2 has been called before. +/// This function is automatically called when the histogram is created +/// if the static function TH1::SetDefaultSumw2 has been called before. /// If flag is false the structure is deleted void TProfile::Sumw2(Bool_t flag) diff --git a/hist/hist/src/TProfile2D.cxx b/hist/hist/src/TProfile2D.cxx index 611d476a76a..d221fc50659 100644 --- a/hist/hist/src/TProfile2D.cxx +++ b/hist/hist/src/TProfile2D.cxx @@ -23,50 +23,51 @@ Bool_t TProfile2D::fgApproximate = kFALSE; ClassImp(TProfile2D) -//______________________________________________________________________________ -// -// Profile2D histograms are used to display the mean -// value of Z and its RMS for each cell in X,Y. -// Profile2D histograms are in many cases an -// elegant replacement of three-dimensional histograms : the inter-relation of three -// measured quantities X, Y and Z can always be visualized by a three-dimensional -// histogram or scatter-plot; its representation on the line-printer is not particularly -// satisfactory, except for sparse data. If Z is an unknown (but single-valued) -// approximate function of X,Y this function is displayed by a profile2D histogram with -// much better precision than by a scatter-plot. -// -// The following formulae show the cumulated contents (capital letters) and the values -// displayed by the printing or plotting routines (small letters) of the elements for cell I, J. -// -// 2 -// H(I,J) = sum Z E(I,J) = sum Z -// l(I,J) = sum l L(I,J) = sum l -// h(I,J) = H(I,J)/L(I,J) s(I,J) = sqrt(E(I,J)/L(I,J)- h(I,J)**2) -// e(I,J) = s(I,J)/sqrt(L(I,J)) -// -// In the special case where s(I,J) is zero (eg, case of 1 entry only in one cell) -// the bin error e(I,J) is computed from the average of the s(I,J) for all cells -// if the static function TProfile2D::Approximate has been called. -// This simple/crude approximation was suggested in order to keep the cell -// during a fit operation. But note that this approximation is not the default behaviour. -// -// Example of a profile2D histogram -//{ -// TCanvas *c1 = new TCanvas("c1","Profile histogram example",200,10,700,500); -// hprof2d = new TProfile2D("hprof2d","Profile of pz versus px and py",40,-4,4,40,-4,4,0,20); -// Float_t px, py, pz; -// for ( Int_t i=0; i<25000; i++) { -// gRandom->Rannor(px,py); -// pz = px*px + py*py; -// hprof2d->Fill(px,py,pz,1); -// } -// hprof2d->Draw(); -//} -// +/** \class TProfile2D + \ingroup Hist + Profile2D histograms are used to display the mean + value of Z and its RMS for each cell in X,Y. + Profile2D histograms are in many cases an + elegant replacement of three-dimensional histograms : the inter-relation of three + measured quantities X, Y and Z can always be visualized by a three-dimensional + histogram or scatter-plot; its representation on the line-printer is not particularly + satisfactory, except for sparse data. If Z is an unknown (but single-valued) + approximate function of X,Y this function is displayed by a profile2D histogram with + much better precision than by a scatter-plot. + + The following formulae show the cumulated contents (capital letters) and the values + displayed by the printing or plotting routines (small letters) of the elements for cell I, J. + + 2 + H(I,J) = sum Z E(I,J) = sum Z + l(I,J) = sum l L(I,J) = sum l + h(I,J) = H(I,J)/L(I,J) s(I,J) = sqrt(E(I,J)/L(I,J)- h(I,J)**2) + e(I,J) = s(I,J)/sqrt(L(I,J)) + + In the special case where s(I,J) is zero (eg, case of 1 entry only in one cell) + the bin error e(I,J) is computed from the average of the s(I,J) for all cells + if the static function TProfile2D::Approximate has been called. + This simple/crude approximation was suggested in order to keep the cell + during a fit operation. But note that this approximation is not the default behaviour. + + Example of a profile2D histogram + ~~~~{.cpp} + { + TCanvas *c1 = new TCanvas("c1","Profile histogram example",200,10,700,500); + hprof2d = new TProfile2D("hprof2d","Profile of pz versus px and py",40,-4,4,40,-4,4,0,20); + Float_t px, py, pz; + for ( Int_t i=0; i<25000; i++) { + gRandom->Rannor(px,py); + pz = px*px + py*py; + hprof2d->Fill(px,py,pz,1); + } + hprof2d->Draw(); + } + ~~~~ +*/ //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Default constructor for Profile2D histograms*-*-*-*-*-*-*-*-* -///*-* ============================================ +/// Default constructor for Profile2D histograms TProfile2D::TProfile2D() : TH2D() { @@ -76,16 +77,14 @@ TProfile2D::TProfile2D() : TH2D() } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Default destructor for Profile2D histograms*-*-*-*-*-*-*-*-* -///*-* =========================================== +/// Default destructor for Profile2D histograms TProfile2D::~TProfile2D() { } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Normal Constructor for Profile histograms*-*-*-*-*-*-*-*-*-* -///*-* ========================================== +/// Normal Constructor for Profile histograms -* /// /// The first eight parameters are similar to TH2D::TH2D. /// All values of z are accepted at filling time. @@ -102,7 +101,7 @@ TProfile2D::~TProfile2D() /// H(I,J) is printed as the cell contents. The errors computed are s(I,J) if CHOPT='S' /// (spread option), or e(I,J) if CHOPT=' ' (error on mean). /// -/// See TProfile2D::BuildOptions for explanation of parameters +/// See TProfile2D::BuildOptions for explanation of parameters /// /// see other constructors below with all possible combinations of /// fix and variable bin size like in TH2D. @@ -143,14 +142,13 @@ TProfile2D::TProfile2D(const char *name,const char *title,Int_t nx,const Double_ //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Constructor for Profile2D histograms with range in z*-*-*-*-*-* -///*-* ==================================================== -/// The first eight parameters are similar to TH2D::TH2D. -/// Only the values of Z between ZMIN and ZMAX will be considered at filling time. -/// zmin and zmax will also be the maximum and minimum values -/// on the z scale when drawing the profile2D. +/// Constructor for Profile2D histograms with range in z +/// The first eight parameters are similar to TH2D::TH2D. +/// Only the values of Z between ZMIN and ZMAX will be considered at filling time. +/// zmin and zmax will also be the maximum and minimum values +/// on the z scale when drawing the profile2D. /// -/// See TProfile2D::BuildOptions for more explanations on errors +/// See TProfile2D::BuildOptions for more explanations on errors /// TProfile2D::TProfile2D(const char *name,const char *title,Int_t nx,Double_t xlow,Double_t xup,Int_t ny, Double_t ylow,Double_t yup,Double_t zlow,Double_t zup,Option_t *option) @@ -162,8 +160,7 @@ TProfile2D::TProfile2D(const char *name,const char *title,Int_t nx,Double_t xlow //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Set Profile2D histogram structure and options*-*-*-*-*-*-*-*-* -///*-* ============================================= +/// Set Profile2D histogram structure and options /// /// zmin: minimum value allowed for z /// zmax: maximum value allowed for z @@ -226,11 +223,9 @@ Bool_t TProfile2D::Add(const TH1 *h1, Double_t c1) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Replace contents of this profile2D by the addition of h1 and h2*-*-* -///*-* =============================================================== +/// Replace contents of this profile2D by the addition of h1 and h2* /// /// this = c1*h1 + c2*h2 -/// Bool_t TProfile2D::Add(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2) { @@ -251,7 +246,7 @@ Bool_t TProfile2D::Add(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2) //////////////////////////////////////////////////////////////////////////////// -/// static function +/// static function /// set the fgApproximate flag. When the flag is true, the function GetBinError /// will approximate the bin error with the average profile error on all bins /// in the following situation only @@ -268,9 +263,9 @@ void TProfile2D::Approximate(Bool_t approx) //////////////////////////////////////////////////////////////////////////////// /// Fill histogram with all entries in the buffer. -/// action = -1 histogram is reset and refilled from the buffer (called by THistPainter::Paint) -/// action = 0 histogram is filled from the buffer -/// action = 1 histogram is filled and buffer is deleted +/// - action = -1 histogram is reset and refilled from the buffer (called by THistPainter::Paint) +/// - action = 0 histogram is filled from the buffer +/// - action = 1 histogram is filled and buffer is deleted /// The buffer is automatically deleted when the number of entries /// in the buffer is greater than the number of entries in the histogram @@ -332,11 +327,11 @@ Int_t TProfile2D::BufferEmpty(Int_t action) //////////////////////////////////////////////////////////////////////////////// /// accumulate arguments in buffer. When buffer is full, empty the buffer -/// fBuffer[0] = number of entries in buffer -/// fBuffer[1] = w of first entry -/// fBuffer[2] = x of first entry -/// fBuffer[3] = y of first entry -/// fBuffer[4] = z of first entry +/// fBuffer[0] = number of entries in buffer +/// fBuffer[1] = w of first entry +/// fBuffer[2] = x of first entry +/// fBuffer[3] = y of first entry +/// fBuffer[4] = z of first entry Int_t TProfile2D::BufferFill(Double_t x, Double_t y, Double_t z, Double_t w) { @@ -364,8 +359,7 @@ Int_t TProfile2D::BufferFill(Double_t x, Double_t y, Double_t z, Double_t w) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*Copy a Profile2D histogram to a new profile2D histogram*-*-*-* -///*-* ======================================================= +/// Copy a Profile2D histogram to a new profile2D histogram* void TProfile2D::Copy(TObject &obj) const { @@ -404,8 +398,7 @@ Bool_t TProfile2D::Divide(TF1 *, Double_t ) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Divide this profile2D by h1*-*-*-*-*-*-*-*-*-*-*-*-* -///*-* =========================== +/// Divide this profile2D by h1 /// /// this = this/h1 /// @@ -489,8 +482,7 @@ Bool_t TProfile2D::Divide(const TH1 *h1) //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Replace contents of this profile2D by the division of h1 by h2*-*-* -///*-* ============================================================== +/// Replace contents of this profile2D by the division of h1 by h2 /// /// this = c1*h1/(c2*h2) /// @@ -591,7 +583,7 @@ Bool_t TProfile2D::Divide(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2 } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Fill a Profile2D histogram (no weights)*-*-*-*-*-*-*-* +/// Fill a Profile2D histogram (no weights) ///*-* ======================================= Int_t TProfile2D::Fill(Double_t x, Double_t y, Double_t z) @@ -745,8 +737,7 @@ Int_t TProfile2D::Fill(const char *namex, Double_t y, Double_t z) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Fill a Profile2D histogram with weights*-*-*-*-*-*-*-* -///*-* ======================================= +/// Fill a Profile2D histogram with weights Int_t TProfile2D::Fill(Double_t x, Double_t y, Double_t z, Double_t w) { @@ -788,8 +779,7 @@ Int_t TProfile2D::Fill(Double_t x, Double_t y, Double_t z, Double_t w) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Return bin content of a Profile2D histogram*-*-*-*-*-*-*-*-* -///*-* =========================================== +/// Return bin content of a Profile2D histogram Double_t TProfile2D::GetBinContent(Int_t bin) const { @@ -802,8 +792,7 @@ Double_t TProfile2D::GetBinContent(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Return bin entries of a Profile2D histogram*-*-*-*-*-*-*-*-* -///*-* =========================================== +/// Return bin entries of a Profile2D histogram Double_t TProfile2D::GetBinEntries(Int_t bin) const { @@ -814,12 +803,12 @@ Double_t TProfile2D::GetBinEntries(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -/// Return bin effective entries for a weighted filled Profile histogram. -/// In case of an unweighted profile, it is equivalent to the number of entries per bin -/// The effective entries is defined as the square of the sum of the weights divided by the -/// sum of the weights square. -/// TProfile::Sumw2() must be called before filling the profile with weights. -/// Only by calling this method the sum of the square of the weights per bin is stored. +/// Return bin effective entries for a weighted filled Profile histogram. +/// In case of an unweighted profile, it is equivalent to the number of entries per bin +/// The effective entries is defined as the square of the sum of the weights divided by the +/// sum of the weights square. +/// TProfile::Sumw2() must be called before filling the profile with weights. +/// Only by calling this method the sum of the square of the weights per bin is stored. /// ///*-* ========================================= @@ -829,7 +818,7 @@ Double_t TProfile2D::GetBinEffectiveEntries(Int_t bin) } //////////////////////////////////////////////////////////////////////////////// -/// *-*-*-*-*-*-*Return bin error of a Profile2D histogram*-*-*-*-*-*-*-*-* +/// -*Return bin error of a Profile2D histogram /// /// Computing errors: A moving field /// ================================= @@ -852,7 +841,7 @@ Double_t TProfile2D::GetBinError(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Return option to compute profile2D errors*-*-*-*-*-*-*-* +/// -*Return option to compute profile2D errors ///*-* ========================================= Option_t *TProfile2D::GetErrorOption() const @@ -1148,7 +1137,7 @@ Bool_t TProfile2D::Multiply(TF1 *, Double_t ) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Multiply this profile2D by h1*-*-*-*-*-*-*-*-*-*-*-* +/// Multiply this profile2D by h1 - ///*-* ============================= /// /// this = this*h1 @@ -1162,7 +1151,7 @@ Bool_t TProfile2D::Multiply(const TH1 *) //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Replace contents of this profile2D by multiplication of h1 by h2*-* +///*-*Replace contents of this profile2D by multiplication of h1 by h2*-* ///*-* ================================================================ /// /// this = (c1*h1)*(c2*h2) @@ -1175,7 +1164,7 @@ Bool_t TProfile2D::Multiply(const TH1 *, const TH1 *, Double_t, Double_t, Option } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Project this profile2D into a 2-D histogram along X,Y*-*-*-*-*-*-* +///*-*Project this profile2D into a 2-D histogram along X,Y -* ///*-* ===================================================== /// /// The projection is always of the type TH2D. @@ -1267,7 +1256,7 @@ TH2D *TProfile2D::ProjectionXY(const char *name, Option_t *option) const } //////////////////////////////////////////////////////////////////////////////// -/// *-*-*-*-*Project a 2-D histogram into a profile histogram along X*-*-*-*-*-* +/// *-*Project a 2-D histogram into a profile histogram along X /// *-* ======================================================== /// /// The projection is made from the channels along the Y axis @@ -1289,7 +1278,7 @@ TProfile *TProfile2D::ProfileX(const char *name, Int_t firstybin, Int_t lastybin } //////////////////////////////////////////////////////////////////////////////// -/// *-*-*-*-*Project a 2-D histogram into a profile histogram along X*-*-*-*-*-* +/// *-*Project a 2-D histogram into a profile histogram along X /// *-* ======================================================== /// /// The projection is made from the channels along the X axis @@ -1404,7 +1393,7 @@ void TProfile2D::PutStats(Double_t *stats) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Reset contents of a Profile2D histogram*-*-*-*-*-*-*-* +/// -*Reset contents of a Profile2D histogram ///*-* ======================================= void TProfile2D::Reset(Option_t *option) @@ -1440,27 +1429,26 @@ void TProfile2D::ExtendAxis(Double_t x, TAxis *axis) } //////////////////////////////////////////////////////////////////////////////// -/// -*-*-*Rebin this histogram grouping nxgroup/nygroup bins along the xaxis/yaxis together*-*-*-*- -/// ================================================================================= -/// if newname is not blank a new profile hnew is created. -/// else the current histogram is modified (default) -/// The parameter nxgroup/nygroup indicate how many bins along the xaxis/yaxis of this -/// have to be merged into one bin of hnew -/// If the original profile has errors stored (via Sumw2), the resulting -/// profile has new errors correctly calculated. +/// Rebin this histogram grouping nxgroup/nygroup bins along the xaxis/yaxis together. +/// if newname is not blank a new profile hnew is created. +/// else the current histogram is modified (default) +/// The parameter nxgroup/nygroup indicate how many bins along the xaxis/yaxis of this +/// have to be merged into one bin of hnew +/// If the original profile has errors stored (via Sumw2), the resulting +/// profile has new errors correctly calculated. /// -/// examples: if hpxpy is an existing TProfile2D profile with 40 x 40 bins -/// hpxpy->Rebin2D(); // merges two bins along the xaxis and yaxis in one -/// // Carefull: previous contents of hpxpy are lost -/// hpxpy->Rebin2D(3,5); // merges 3 bins along the xaxis and 5 bins along the yaxis in one -/// // Carefull: previous contents of hpxpy are lost -/// hpxpy->RebinX(5); //merges five bins along the xaxis in one in hpxpy -/// TProfile2D *hnew = hpxpy->RebinY(5,"hnew"); // creates a new profile hnew -/// // merging 5 bins of hpxpy along the yaxis in one bin +/// examples: if hpxpy is an existing TProfile2D profile with 40 x 40 bins +/// hpxpy->Rebin2D(); // merges two bins along the xaxis and yaxis in one +/// // Carefull: previous contents of hpxpy are lost +/// hpxpy->Rebin2D(3,5); // merges 3 bins along the xaxis and 5 bins along the yaxis in one +/// // Carefull: previous contents of hpxpy are lost +/// hpxpy->RebinX(5); //merges five bins along the xaxis in one in hpxpy +/// TProfile2D *hnew = hpxpy->RebinY(5,"hnew"); // creates a new profile hnew +/// // merging 5 bins of hpxpy along the yaxis in one bin /// -/// NOTE : If nxgroup/nygroup is not an exact divider of the number of bins, -/// along the xaxis/yaxis the top limit(s) of the rebinned profile -/// is changed to the upper edge of the xbin=newxbins*nxgroup resp. +/// NOTE : If nxgroup/nygroup is not an exact divider of the number of bins, +/// along the xaxis/yaxis the top limit(s) of the rebinned profile +/// is changed to the upper edge of the xbin=newxbins*nxgroup resp. /// ybin=newybins*nygroup and the remaining bins are added to /// the overflow bin. /// Statistics will be recomputed from the new bin contents. @@ -1818,13 +1806,11 @@ void TProfile2D::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/) } //////////////////////////////////////////////////////////////////////////////// -/// *-*-*-*-*Multiply this profile2D by a constant c1*-*-*-*-*-*-*-*-* -/// *-* ======================================== +/// Multiply this profile2D by a constant c1 /// /// this = c1*this /// /// This function uses the services of TProfile2D::Add -/// void TProfile2D::Scale(Double_t c1, Option_t * option) { @@ -1832,8 +1818,7 @@ void TProfile2D::Scale(Double_t c1, Option_t * option) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*Set the number of entries in bin*-*-*-*-*-*-*-*-*-*-*-* -///*-* ================================ +/// Set the number of entries in bin - void TProfile2D::SetBinEntries(Int_t bin, Double_t w) { @@ -1841,8 +1826,7 @@ void TProfile2D::SetBinEntries(Int_t bin, Double_t w) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*Redefine x and y axis parameters*-*-*-*-*-*-*-*-*-*-*-* -///*-* =========================== +/// Redefine x and y axis parameters - void TProfile2D::SetBins(Int_t nx, Double_t xmin, Double_t xmax, Int_t ny, Double_t ymin, Double_t ymax) { @@ -1852,8 +1836,7 @@ void TProfile2D::SetBins(Int_t nx, Double_t xmin, Double_t xmax, Int_t ny, Doubl } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*Redefine x and y axis parameters for variable bin sizes -*-*-*-*-*-*-* -///*-* =========================== +/// Redefine x and y axis parameters for variable bin sizes - -* void TProfile2D::SetBins(Int_t nx, const Double_t *xbins, Int_t ny, const Double_t *ybins) { @@ -1893,24 +1876,19 @@ void TProfile2D::SetBuffer(Int_t buffersize, Option_t *) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Set option to compute profile2D errors*-*-*-*-*-*-*-* -///*-* ======================================= -/// The computation of the bin errors is based on the parameter option: -/// option: -/// ' ' (Default) The bin errors are the standard error on the mean of the bin profiled values (Z), +/// Set option to compute profile2D errors +/// The computation of the bin errors is based on the parameter option: +/// - ' ' (Default) The bin errors are the standard error on the mean of the bin profiled values (Z), /// i.e. the standard error of the bin contents. /// Note that if TProfile::Approximate() is called, an approximation is used when /// the spread in Z is 0 and the number of bin entries is > 0 -/// -/// 's' The bin errors are the standard deviations of the Z bin values +/// - 's' The bin errors are the standard deviations of the Z bin values /// Note that if TProfile::Approximate() is called, an approximation is used when /// the spread in Z is 0 and the number of bin entries is > 0 -/// -/// 'i' Errors are as in default case (standard errors of the bin contents) +/// - 'i' Errors are as in default case (standard errors of the bin contents) /// The only difference is for the case when the spread in Z is zero. /// In this case for N > 0 the error is 1./SQRT(12.*N) -/// -/// 'g' Errors are 1./SQRT(W) for W not equal to 0 and 0 for W = 0. +/// - 'g' Errors are 1./SQRT(W) for W not equal to 0 and 0 for W = 0. /// W is the sum in the bin of the weights of the profile. /// This option is for combining measurements z +/- dz, /// and the profile is filled with values y and weights z = 1/dz**2 @@ -1957,10 +1935,9 @@ void TProfile2D::Streamer(TBuffer &R__b) } //////////////////////////////////////////////////////////////////////////////// -/// Create/Delete structure to store sum of squares of weights per bin *-*-*-*-*-*-*-* -/// This is needed to compute the correct statistical quantities -/// of a profile filled with weights -/// +/// Create/Delete structure to store sum of squares of weights per bin +/// This is needed to compute the correct statistical quantities +/// of a profile filled with weights /// /// This function is automatically called when the histogram is created /// if the static function TH1::SetDefaultSumw2 has been called before. diff --git a/hist/hist/src/TProfile3D.cxx b/hist/hist/src/TProfile3D.cxx index 87b96192ace..d3e9840eb44 100644 --- a/hist/hist/src/TProfile3D.cxx +++ b/hist/hist/src/TProfile3D.cxx @@ -25,53 +25,54 @@ Bool_t TProfile3D::fgApproximate = kFALSE; ClassImp(TProfile3D) -//______________________________________________________________________________ -// -// Profile3D histograms are used to display the mean -// value of T and its RMS for each cell in X,Y,Z. -// Profile3D histograms are in many cases an -// The inter-relation of three measured quantities X, Y, Z and T can always -// be visualized by a four-dimensional histogram or scatter-plot; -// its representation on the line-printer is not particularly -// satisfactory, except for sparse data. If T is an unknown (but single-valued) -// approximate function of X,Y,Z this function is displayed by a profile3D histogram with -// much better precision than by a scatter-plot. -// -// The following formulae show the cumulated contents (capital letters) and the values -// displayed by the printing or plotting routines (small letters) of the elements for cell I, J. -// -// 2 -// H(I,J,K) = sum T E(I,J,K) = sum T -// l(I,J,K) = sum l L(I,J,K) = sum l -// h(I,J,K) = H(I,J,K)/L(I,J,K) s(I,J,K) = sqrt(E(I,J,K)/L(I,J,K)- h(I,J,K)**2) -// e(I,J,K) = s(I,J,K)/sqrt(L(I,J,K)) -// -// In the special case where s(I,J,K) is zero (eg, case of 1 entry only in one cell) -// e(I,J,K) is computed from the average of the s(I,J,K) for all cells, -// if the static function TProfile3D::Approximate has been called. -// This simple/crude approximation was suggested in order to keep the cell -// during a fit operation. But note that this approximation is not the default behaviour. -// -// Example of a profile3D histogram -//{ -// TCanvas *c1 = new TCanvas("c1","Profile histogram example",200,10,700,500); -// hprof3d = new TProfile3D("hprof3d","Profile of pt versus px, py and pz",40,-4,4,40,-4,4,40,0,20); -// Double_t px, py, pz, pt; -// TRandom3 r(0); -// for ( Int_t i=0; i<25000; i++) { -// r.Rannor(px,py); -// pz = px*px + py*py; -// pt = r.Landau(0,1); -// hprof3d->Fill(px,py,pz,pt,1); -// } -// hprof3d->Draw(); -//} -// -// NOTE: A TProfile3D is drawn as it was a simple TH3 +/** \class TProfile3D + \ingroup Hist + Profile3D histograms are used to display the mean + value of T and its RMS for each cell in X,Y,Z. + Profile3D histograms are in many cases an + The inter-relation of three measured quantities X, Y, Z and T can always + be visualized by a four-dimensional histogram or scatter-plot; + its representation on the line-printer is not particularly + satisfactory, except for sparse data. If T is an unknown (but single-valued) + approximate function of X,Y,Z this function is displayed by a profile3D histogram with + much better precision than by a scatter-plot. + + The following formulae show the cumulated contents (capital letters) and the values + displayed by the printing or plotting routines (small letters) of the elements for cell I, J. + + 2 + H(I,J,K) = sum T E(I,J,K) = sum T + l(I,J,K) = sum l L(I,J,K) = sum l + h(I,J,K) = H(I,J,K)/L(I,J,K) s(I,J,K) = sqrt(E(I,J,K)/L(I,J,K)- h(I,J,K)**2) + e(I,J,K) = s(I,J,K)/sqrt(L(I,J,K)) + + In the special case where s(I,J,K) is zero (eg, case of 1 entry only in one cell) + e(I,J,K) is computed from the average of the s(I,J,K) for all cells, + if the static function TProfile3D::Approximate has been called. + This simple/crude approximation was suggested in order to keep the cell + during a fit operation. But note that this approximation is not the default behaviour. + + Example of a profile3D histogram +~~~~{.cpp} +{ + TCanvas *c1 = new TCanvas("c1","Profile histogram example",200,10,700,500); + hprof3d = new TProfile3D("hprof3d","Profile of pt versus px, py and pz",40,-4,4,40,-4,4,40,0,20); + Double_t px, py, pz, pt; + TRandom3 r(0); + for ( Int_t i=0; i<25000; i++) { + r.Rannor(px,py); + pz = px*px + py*py; + pt = r.Landau(0,1); + hprof3d->Fill(px,py,pz,pt,1); + } + hprof3d->Draw(); +} +~~~~ + NOTE: A TProfile3D is drawn as it was a simple TH3 +*/ //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Default constructor for Profile3D histograms*-*-*-*-*-*-*-*-* -///*-* ============================================ +/// Default constructor for Profile3D histograms TProfile3D::TProfile3D() : TH3D() { @@ -81,19 +82,17 @@ TProfile3D::TProfile3D() : TH3D() } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Default destructor for Profile3D histograms*-*-*-*-*-*-*-*-* -///*-* =========================================== +/// Default destructor for Profile3D histograms TProfile3D::~TProfile3D() { } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*Normal Constructor for Profile histograms*-*-*-*-*-*-*-*-*-* -///*-* ========================================== +/// Normal Constructor for Profile histograms-* /// -/// The first eleven parameters are similar to TH3D::TH3D. -/// All values of t are accepted at filling time. +/// The first eleven parameters are similar to TH3D::TH3D. +/// All values of t are accepted at filling time. /// To fill a profile3D histogram, one must use TProfile3D::Fill function. /// /// Note that when filling the profile histogram the function Fill @@ -107,10 +106,10 @@ TProfile3D::~TProfile3D() /// H(I,J,K) is printed as the cell contents. The errors computed are s(I,J,K) if CHOPT='S' /// (spread option), or e(I,J,K) if CHOPT=' ' (error on mean). /// -/// See TProfile3D::BuildOptions for explanation of parameters +/// See TProfile3D::BuildOptions for explanation of parameters /// -/// see other constructors below with all possible combinations of -/// fix and variable bin size like in TH3D. +/// see other constructors below with all possible combinations of +/// fix and variable bin size like in TH3D. TProfile3D::TProfile3D(const char *name,const char *title,Int_t nx,Double_t xlow,Double_t xup,Int_t ny,Double_t ylow,Double_t yup,Int_t nz, Double_t zlow,Double_t zup,Option_t *option) : TH3D(name,title,nx,xlow,xup,ny,ylow,yup,nz,zlow,zup) @@ -129,8 +128,7 @@ TProfile3D::TProfile3D(const char *name,const char *title,Int_t nx,const Double_ } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Set Profile3D histogram structure and options*-*-*-*-*-*-*-*-* -///*-* ============================================= +/// Set Profile3D histogram structure and options /// /// tmin: minimum value allowed for t /// tmax: maximum value allowed for t @@ -155,7 +153,7 @@ void TProfile3D::BuildOptions(Double_t tmin, Double_t tmax, Option_t *option) } //////////////////////////////////////////////////////////////////////////////// -///copy constructor +/// copy constructor TProfile3D::TProfile3D(const TProfile3D &profile) : TH3D() { @@ -191,11 +189,9 @@ Bool_t TProfile3D::Add(const TH1 *h1, Double_t c1) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Replace contents of this profile3D by the addition of h1 and h2*-*-* -///*-* =============================================================== +/// Replace contents of this profile3D by the addition of h1 and h2 /// /// this = c1*h1 + c2*h2 -/// Bool_t TProfile3D::Add(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2) { @@ -217,7 +213,6 @@ Bool_t TProfile3D::Add(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2) //////////////////////////////////////////////////////////////////////////////// -/// static function /// set the fgApproximate flag. When the flag is true, the function GetBinError /// will approximate the bin error with the average profile error on all bins /// in the following situation only @@ -234,9 +229,9 @@ void TProfile3D::Approximate(Bool_t approx) //////////////////////////////////////////////////////////////////////////////// /// Fill histogram with all entries in the buffer. -/// action = -1 histogram is reset and refilled from the buffer (called by THistPainter::Paint) -/// action = 0 histogram is filled from the buffer -/// action = 1 histogram is filled and buffer is deleted +/// - action = -1 histogram is reset and refilled from the buffer (called by THistPainter::Paint) +/// - action = 0 histogram is filled from the buffer +/// - action = 1 histogram is filled and buffer is deleted /// The buffer is automatically deleted when the number of entries /// in the buffer is greater than the number of entries in the histogram @@ -339,8 +334,7 @@ Int_t TProfile3D::BufferFill(Double_t x, Double_t y, Double_t z, Double_t t, Dou } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*Copy a Profile3D histogram to a new profile2D histogram*-*-*-* -///*-* ======================================================= +/// Copy a Profile3D histogram to a new profile2D histogram void TProfile3D::Copy(TObject &obj) const { @@ -377,8 +371,7 @@ Bool_t TProfile3D::Divide(TF1 *, Double_t ) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Divide this profile2D by h1*-*-*-*-*-*-*-*-*-*-*-*-* -///*-* =========================== +/// Divide this profile2D by h1 /// /// this = this/h1 /// @@ -473,8 +466,7 @@ Bool_t TProfile3D::Divide(const TH1 *h1) //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Replace contents of this profile2D by the division of h1 by h2*-*-* -///*-* ============================================================== +/// Replace contents of this profile2D by the division of h1 by h2 /// /// this = c1*h1/(c2*h2) /// @@ -584,8 +576,7 @@ Bool_t TProfile3D::Divide(const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2 } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Fill a Profile3D histogram (no weights)*-*-*-*-*-*-*-* -///*-* ======================================= +/// Fill a Profile3D histogram (no weights) Int_t TProfile3D::Fill(Double_t x, Double_t y, Double_t z, Double_t t) { @@ -634,8 +625,7 @@ Int_t TProfile3D::Fill(Double_t x, Double_t y, Double_t z, Double_t t) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Fill a Profile3D histogram with weights*-*-*-*-*-*-*-* -///*-* ======================================= +/// Fill a Profile3D histogram with weights Int_t TProfile3D::Fill(Double_t x, Double_t y, Double_t z, Double_t t, Double_t w) { @@ -685,8 +675,7 @@ Int_t TProfile3D::Fill(Double_t x, Double_t y, Double_t z, Double_t t, Double_t } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Return bin content of a Profile3D histogram*-*-*-*-*-*-*-*-* -///*-* =========================================== +/// Return bin content of a Profile3D histogram Double_t TProfile3D::GetBinContent(Int_t bin) const { @@ -699,8 +688,7 @@ Double_t TProfile3D::GetBinContent(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*Return bin entries of a Profile3D histogram*-*-*-*-*-*-*-*-* -///*-* =========================================== +/// Return bin entries of a Profile3D histogram Double_t TProfile3D::GetBinEntries(Int_t bin) const { @@ -711,12 +699,12 @@ Double_t TProfile3D::GetBinEntries(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -/// Return bin effective entries for a weighted filled Profile histogram. -/// In case of an unweighted profile, it is equivalent to the number of entries per bin -/// The effective entries is defined as the square of the sum of the weights divided by the -/// sum of the weights square. -/// TProfile::Sumw2() must be called before filling the profile with weights. -/// Only by calling this method the sum of the square of the weights per bin is stored. +/// Return bin effective entries for a weighted filled Profile histogram. +/// In case of an unweighted profile, it is equivalent to the number of entries per bin +/// The effective entries is defined as the square of the sum of the weights divided by the +/// sum of the weights square. +/// TProfile::Sumw2() must be called before filling the profile with weights. +/// Only by calling this method the sum of the square of the weights per bin is stored. /// ///*-* ========================================= @@ -726,7 +714,7 @@ Double_t TProfile3D::GetBinEffectiveEntries(Int_t bin) } //////////////////////////////////////////////////////////////////////////////// -/// *-*-*-*-*-*-*Return bin error of a Profile3D histogram*-*-*-*-*-*-*-*-* +/// -*Return bin error of a Profile3D histogram /// /// Computing errors: A moving field /// ================================= @@ -749,7 +737,7 @@ Double_t TProfile3D::GetBinError(Int_t bin) const } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Return option to compute profile2D errors*-*-*-*-*-*-*-* +///-*Return option to compute profile2D errors ///*-* ========================================= Option_t *TProfile3D::GetErrorOption() const @@ -1029,7 +1017,7 @@ Bool_t TProfile3D::Multiply(TF1 *, Double_t ) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*-*Multiply this profile2D by h1*-*-*-*-*-*-*-*-*-*-*-* +///-*-*Multiply this profile2D by h1 - ///*-* ============================= /// /// this = this*h1 @@ -1056,7 +1044,7 @@ Bool_t TProfile3D::Multiply(const TH1 *, const TH1 *, Double_t, Double_t, Option } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*Project this profile3D into a 3-D histogram along X,Y,Z*-*-*-*-*-*-* +///*-*-*-*-*Project this profile3D into a 3-D histogram along X,Y,Z -* ///*-* ===================================================== /// /// The projection is always of the type TH3D. @@ -1165,7 +1153,7 @@ TH3D *TProfile3D::ProjectionXYZ(const char *name, Option_t *option) const /// option = "zx" return the z versus x projection into a TProfile2D histogram /// option = "yz" return the y versus z projection into a TProfile2D histogram /// option = "zy" return the z versus y projection into a TProfile2D histogram -/// NB: the notation "a vs b" means "a" vertical and "b" horizontalalong X*-*-*-*-*-* +/// NB: the notation "a vs b" means "a" vertical and "b" horizontalalong X /// /// The resulting profile contains the combination of all the considered bins along X /// By default, all bins are included considering also underflow/overflows @@ -1288,7 +1276,7 @@ void TProfile3D::PutStats(Double_t *stats) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Reset contents of a Profile3D histogram*-*-*-*-*-*-*-* +///-*Reset contents of a Profile3D histogram ///*-* ======================================= void TProfile3D::Reset(Option_t *option) @@ -1381,13 +1369,11 @@ void TProfile3D::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/) } //////////////////////////////////////////////////////////////////////////////// -/// *-*-*-*-*Multiply this profile2D by a constant c1*-*-*-*-*-*-*-*-* -/// *-* ======================================== +/// Multiply this profile2D by a constant c1 /// /// this = c1*this /// /// This function uses the services of TProfile3D::Add -/// void TProfile3D::Scale(Double_t c1, Option_t *option) { @@ -1395,8 +1381,7 @@ void TProfile3D::Scale(Double_t c1, Option_t *option) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*Set the number of entries in bin*-*-*-*-*-*-*-*-*-*-*-* -///*-* ================================ +///Set the number of entries in bin void TProfile3D::SetBinEntries(Int_t bin, Double_t w) { @@ -1404,8 +1389,7 @@ void TProfile3D::SetBinEntries(Int_t bin, Double_t w) } //////////////////////////////////////////////////////////////////////////////// -/// -*-*-*-*-*-*-*Redefine x, y and z axis parameters*-*-*-*-*-*-*-*-*-*-*-* -///*-* =========================== +/// Redefine x, y and z axis parameters void TProfile3D::SetBins(Int_t nx, Double_t xmin, Double_t xmax, Int_t ny, Double_t ymin, Double_t ymax, Int_t nz, Double_t zmin, Double_t zmax) { @@ -1415,8 +1399,7 @@ void TProfile3D::SetBins(Int_t nx, Double_t xmin, Double_t xmax, Int_t ny, Doubl } //////////////////////////////////////////////////////////////////////////////// -/// -*-*-*-*-*-*-*Redefine x, y and z axis parameters with variable bin sizes *-*-*-*-*-*-*-*-* -/// ============================================================ +/// Redefine x, y and z axis parameters with variable bin sizes void TProfile3D::SetBins(Int_t nx, const Double_t *xBins, Int_t ny, const Double_t *yBins, Int_t nz, const Double_t *zBins) { @@ -1456,25 +1439,20 @@ void TProfile3D::SetBuffer(Int_t buffersize, Option_t *) } //////////////////////////////////////////////////////////////////////////////// -///*-*-*-*-*-*-*-*-*-*Set option to compute profile3D errors*-*-*-*-*-*-*-* -///*-* ======================================= +/// Set option to compute profile3D errors /// /// The computation of the bin errors is based on the parameter option: -/// option: -/// ' ' (Default) The bin errors are the standard error on the mean of the bin profiled values (T), +/// - ' ' (Default) The bin errors are the standard error on the mean of the bin profiled values (T), /// i.e. the standard error of the bin contents. /// Note that if TProfile3D::Approximate() is called, an approximation is used when /// the spread in T is 0 and the number of bin entries is > 0 -/// -/// 's' The bin errors are the standard deviations of the T bin values +/// - 's' The bin errors are the standard deviations of the T bin values /// Note that if TProfile3D::Approximate() is called, an approximation is used when /// the spread in T is 0 and the number of bin entries is > 0 -/// -/// 'i' Errors are as in default case (standard errors of the bin contents) +/// - 'i' Errors are as in default case (standard errors of the bin contents) /// The only difference is for the case when the spread in T is zero. /// In this case for N > 0 the error is 1./SQRT(12.*N) -/// -/// 'g' Errors are 1./SQRT(W) for W not equal to 0 and 0 for W = 0. +/// - 'g' Errors are 1./SQRT(W) for W not equal to 0 and 0 for W = 0. /// W is the sum in the bin of the weights of the profile. /// This option is for combining measurements t +/- dt, /// and the profile is filled with values t and weights w = 1/dt**2 @@ -1487,9 +1465,9 @@ void TProfile3D::SetErrorOption(Option_t *option) } //////////////////////////////////////////////////////////////////////////////// -/// Create/Delete structure to store sum of squares of weights per bin *-*-*-*-*-*-*-* -/// This is needed to compute the correct statistical quantities -/// of a profile filled with weights +/// Create/Delete structure to store sum of squares of weights per bin +/// This is needed to compute the correct statistical quantities +/// of a profile filled with weights /// /// /// This function is automatically called when the histogram is created diff --git a/hist/hist/src/TSVDUnfold.cxx b/hist/hist/src/TSVDUnfold.cxx index 4d26d78ddc5..ceda9ccd795 100644 --- a/hist/hist/src/TSVDUnfold.cxx +++ b/hist/hist/src/TSVDUnfold.cxx @@ -25,8 +25,9 @@ //////////////////////////////////////////////////////////////////////////////// -/* Begin_Html - <center><h2>SVD Approach to Data Unfolding</h2></center> +/** \class TSVDUnfold + \ingroup Hist + SVD Approach to Data Unfolding <p> Reference: <a href="http://arXiv.org/abs/hep-ph/9509307">Nucl. Instrum. Meth. A372, 469 (1996) [hep-ph/9509307]</a> <p> @@ -56,9 +57,7 @@ Covariance matrices on the measured spectrum (for either the total uncertainties or individual sources of uncertainties) can be propagated to covariance matrices using the <tt>GetUnfoldCovMatrix</tt> method, which uses pseudo experiments for the propagation. In addition, <tt>GetAdetCovMatrix</tt> allows for the propagation of the statistical uncertainties on the response matrix using pseudo experiments. The covariance matrix corresponding to <tt>Bcov</tt> is also computed as described in <a href="http://arXiv.org/abs/hep-ph/9509307">Nucl. Instrum. Meth. A372, 469 (1996) [hep-ph/9509307]</a> and can be obtained from <tt>tsvdunf->GetXtau()</tt> and its (regularisation independent) inverse from <tt>tsvdunf->GetXinv()</tt>. The distribution of singular values can be retrieved using <tt>tsvdunf->GetSV()</tt>. <p> See also the tutorial for a toy example. - End_Html */ -//_______________________________________________________________________ - +*/ #include <iostream> diff --git a/hist/hist/src/TSpline.cxx b/hist/hist/src/TSpline.cxx index e896b3494b5..dd9d2e57e81 100644 --- a/hist/hist/src/TSpline.cxx +++ b/hist/hist/src/TSpline.cxx @@ -9,14 +9,10 @@ * For the list of contributors see $ROOTSYS/README/CREDITS. * *************************************************************************/ -////////////////////////////////////////////////////////////////////////// -// // -// TSpline // -// // -// Base class for spline implementation containing the Draw/Paint // -// methods // -// // -////////////////////////////////////////////////////////////////////////// +/** \class TSpline + \ingroup Hist + Base class for spline implementation containing the Draw/Paint methods // +*/ #include "TROOT.h" #include "TSpline.h" diff --git a/hist/hist/src/TUnfold.cxx b/hist/hist/src/TUnfold.cxx index 78c19d0a84d..0ae957c054f 100644 --- a/hist/hist/src/TUnfold.cxx +++ b/hist/hist/src/TUnfold.cxx @@ -25,242 +25,223 @@ // Version 1, added ScanLcurve() method // Version 0, stable version of basic unfolding algorithm -//////////////////////////////////////////////////////////////////////// -// -// TUnfold is used to decompose a measurement y into several sources x -// given the measurement uncertainties and a matrix of migrations A -// -// ************************************************************* -// * For most applications, it is better to use TUnfoldDensity * -// * instead of using TUnfoldSys or TUnfold * -// ************************************************************* -// -// If you use this software, please consider the following citation -// S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201] -// -// More documentation and updates are available on -// http://www.desy.de/~sschmitt -// -// -// A short summary of the algorithm is given in the following: -// -// the "best" x matching the measurement y within errors -// is determined by minimizing the function -// L1+L2+L3 -// -// where -// L1 = (y-Ax)# Vyy^-1 (y-Ax) -// L2 = tau^2 (L(x-x0))# L(x-x0) -// L3 = lambda sum_i(y_i -(Ax)_i) -// -// [the notation # means that the matrix is transposed ] -// -// The term L1 is familiar from a least-square minimisation -// The term L2 defines the regularisation (smootheness condition on x), -// where the parameter tau^2 gives the strength of teh regularisation -// The term L3 is an optional area constraint with Lagrangian parameter -// lambda, ensuring that that the normalisation of x is consistent with the -// normalisation of y -// -// The method can be applied to a very large number of problems, -// where the measured distribution y is a linear superposition -// of several Monte Carlo shapes -// -// Input from measurement: -// ----------------------- -// y: vector of measured quantities (dimension ny) -// Vyy: covariance matrix for y (dimension ny x ny) -// in many cases V is diagonal and calculated from the errors of y -// -// From simulation: -// ---------------- -// A: migration matrix (dimension ny x nx) -// -// Result -// ------ -// x: unknown underlying distribution (dimension nx) -// The error matrix of x, V_xx, is also determined -// -// Regularisation -// -------------- -// tau: parameter, defining the regularisation strength -// L: matrix of regularisation conditions (dimension nl x nx) -// depends on the structure of the input data -// x0: bias distribution, from simulation -// -// Preservation of the area -// ------------------------ -// lambda: lagrangian multiplier -// y_i: one component of the vector y -// (Ax)_i: one component of the vector Ax -// -// -// Determination of the unfolding result x: -// -// (a) not constrained: minimisation is performed as a function of x -// for fixed lambda=0 -// or -// (b) constrained: stationary point is found as a function of x and lambda -// -// The constraint can be useful to reduce biases on the result x -// in cases where the vector y follows non-Gaussian probability densities -// (example: Poisson statistics at counting experiments in particle physics) -// -// Some random examples: -// ====================== -// (1) measure a cross-section as a function of, say, E_T(detector) -// and unfold it to obtain the underlying distribution E_T(generator) -// (2) measure a lifetime distribution and unfold the contributions from -// different flavours -// (3) measure the transverse mass and decay angle -// and unfold for the true mass distribution plus background -// -// Documentation -// ============= -// Some technical documentation is available here: -// http://www.desy.de/~sschmitt -// -// Note: -// For most applications it is better to use the derived class -// TUnfoldSys or even better TUnfoldDensity -// -// TUnfoldSys extends the functionality of TUnfold -// such that systematic errors are propagated to teh result -// and that the unfolding can be done with proper background -// subtraction -// -// TUnfoldDensity extends further the functionality of TUnfoldSys -// complex binning schemes are supported -// The binning of input histograms is handeld consistently: -// (1) the regularisation may be done by density, -// i.e respecting the bin widths -// (2) methods are provided which preserve the proper binning -// of the result histograms -// -// Implementation -// ============== -// The result of the unfolding is calculated as follows: -// -// Lsquared = L#L regularisation conditions squared -// -// epsilon_j = sum_i A_ij vector of efficiencies -// -// E^-1 = ((A# Vyy^-1 A)+tau^2 Lsquared) -// -// x = E (A# Vyy^-1 y + tau^2 Lsquared x0 +lambda/2 * epsilon) is the result -// -// The derivatives -// dx_k/dy_i -// dx_k/dA_ij -// dx_k/d(tau^2) -// are calculated for further usage. -// -// The covariance matrix V_xx is calculated as: -// Vxx_ij = sum_kl dx_i/dy_k Vyy_kl dx_j/dy_l -// -// Warning: -// ======== -// The algorithm is based on "standard" matrix inversion, with the -// known limitations in numerical accuracy and computing cost for -// matrices with large dimensions. -// -// Thus the algorithm should not used for large dimensions of x and y -// nx should not be much larger than 200 -// ny should not be much larger than 1000 -// -// -// Proper choice of tau -// ==================== -// One of the difficult questions is about the choice of tau. -// The method implemented in TUnfold is the L-curve method: -// a two-dimensional curve is plotted -// x-axis: log10(chisquare) -// y-axis: log10(regularisation condition) -// In many cases this curve has an L-shape. The best choice of tau is in the -// kink of the L -// -// Within TUnfold a simple version of the L-curve analysis is available. -// It tests a given number of points in a predefined tau-range and searches -// for the maximum of the curvature in the L-curve (kink position). -// if no tau range is given, the range of the scan is determined automatically -// -// A nice overview of the L-curve method is given in: -// The L-curve and Its Use in the Numerical Treatment of Inverse Problems -// (2000) by P. C. Hansen, in Computational Inverse Problems in -// Electrocardiology, ed. P. Johnston, -// Advances in Computational Bioengineering -// http://www.imm.dtu.dk/~pch/TR/Lcurve.ps -// -// Alternative Regularisation conditions -// ===================================== -// Regularisation is needed for most unfolding problems, in order to avoid -// large oscillations and large correlations on the output bins. -// It means that some extra conditions are applied on the output bins -// -// Within TUnfold these conditions are posed on the difference (x-x0), where -// x: unfolding output -// x0: the bias distribution, by default calculated from -// the input matrix A. There is a method SetBias() to change the -// bias distribution. -// The 3rd argument to DoUnfold() is a scale factor applied to the bias -// bias_default[j] = sum_i A[i][j] -// x0[j] = scaleBias*bias[j] -// The scale factor can be used to -// (a) completely suppress the bias by setting it to zero -// (b) compensate differences in the normalisation between data -// and Monte Carlo -// -// If the regularisation is strong, i.e. large parameter tau, -// then the distribution x or its derivatives will look like the bias -// distribution. If the parameter tau is small, the distribution x is -// independent of the bias. -// -// Three basic types of regularisation are implemented in TUnfold -// -// condition regularisation -// ------------------------------------------------------ -// kRegModeNone none -// kRegModeSize minimize the size of (x-x0) -// kRegModeDerivative minimize the 1st derivative of (x-x0) -// kRegModeCurvature minimize the 2nd derivative of (x-x0) -// -// kRegModeSize is the regularisation scheme which often is found in -// literature. The second derivative is often named curvature. -// Sometimes the bias is not discussed, equivalent to a bias scale factor -// of zero. -// -// The regularisation schemes kRegModeDerivative and -// kRegModeCurvature have the nice feature that they create correlations -// between x-bins, whereas the non-regularized unfolding tends to create -// negative correlations between bins. For these regularisation schemes the -// parameter tau could be tuned such that the correlations are smallest, -// as an alternative to the L-curve method. -// -// If kRegModeSize is chosen or if x is a smooth function through all bins, -// the regularisation condition can be set on all bins together by giving -// the appropriate argument in the constructor (see examples above). -// -// If x is composed of independent groups of bins (for example, -// signal and background binning in two variables), it may be necessary to -// set regularisation conditions for the individual groups of bins. -// In this case, give kRegModeNone in the constructor and specify -// the bin grouping with calls to -// RegularizeBins() specify a 1-dimensional group of bins -// RegularizeBins2D() specify a 2-dimensional group of bins -// -// Note, the class TUnfoldDensity provides an automatic setup of complex -// regularisation schemes -// -// For ultimate flexibility, the regularisation condition can be set on each -// bin individually -// -> give kRegModeNone in the constructor and use -// RegularizeSize() regularize one bin -// RegularizeDerivative() regularize the slope given by two bins -// RegularizeCurvature() regularize the curvature given by three bins -// AddRegularisationCondition() -// define an arbitrary regulatisation condition -// -/////////////////////////////////////////////////////////////////////////// +/** \class TUnfold + \ingroup Hist + TUnfold is used to decompose a measurement y into several sources x + given the measurement uncertainties and a matrix of migrations A + + **For most applications, it is better to use TUnfoldDensity + instead of using TUnfoldSys or TUnfold** + + If you use this software, please consider the following citation + S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201] + + More documentation and updates are available on + http://www.desy.de/~sschmitt + + A short summary of the algorithm is given in the following: + the "best" x matching the measurement y within errors + is determined by minimizing the function + L1+L2+L3 + + where + L1 = (y-Ax)# Vyy^-1 (y-Ax) + L2 = tau^2 (L(x-x0))# L(x-x0) + L3 = lambda sum_i(y_i -(Ax)_i) + + [the notation # means that the matrix is transposed ] + + The term L1 is familiar from a least-square minimisation + The term L2 defines the regularisation (smootheness condition on x), + where the parameter tau^2 gives the strength of teh regularisation + The term L3 is an optional area constraint with Lagrangian parameter + lambda, ensuring that that the normalisation of x is consistent with the + normalisation of y + + The method can be applied to a very large number of problems, + where the measured distribution y is a linear superposition + of several Monte Carlo shapes + + ## Input from measurement: + y: vector of measured quantities (dimension ny) + Vyy: covariance matrix for y (dimension ny x ny) + in many cases V is diagonal and calculated from the errors of y + + ## From simulation: + A: migration matrix (dimension ny x nx) + + ## Result + x: unknown underlying distribution (dimension nx) + The error matrix of x, V_xx, is also determined + + ## Regularisation + tau: parameter, defining the regularisation strength + L: matrix of regularisation conditions (dimension nl x nx) + depends on the structure of the input data + x0: bias distribution, from simulation + + ## Preservation of the area + lambda: lagrangian multiplier + y_i: one component of the vector y + (Ax)_i: one component of the vector Ax + + + Determination of the unfolding result x: + - (a) not constrained: minimisation is performed as a function of x + for fixed lambda=0 + - (b) constrained: stationary point is found as a function of x and lambda + + The constraint can be useful to reduce biases on the result x + in cases where the vector y follows non-Gaussian probability densities + (example: Poisson statistics at counting experiments in particle physics) + + Some random examples: + 1. measure a cross-section as a function of, say, E_T(detector) + and unfold it to obtain the underlying distribution E_T(generator) + 2. measure a lifetime distribution and unfold the contributions from + different flavours + 3. measure the transverse mass and decay angle + and unfold for the true mass distribution plus background + + ## Documentation + Some technical documentation is available here: + http://www.desy.de/~sschmitt + + Note: + + For most applications it is better to use the derived class + TUnfoldSys or even better TUnfoldDensity + + TUnfoldSys extends the functionality of TUnfold + such that systematic errors are propagated to teh result + and that the unfolding can be done with proper background + subtraction + + TUnfoldDensity extends further the functionality of TUnfoldSys + complex binning schemes are supported + The binning of input histograms is handeld consistently: + (1) the regularisation may be done by density, + i.e respecting the bin widths + (2) methods are provided which preserve the proper binning + of the result histograms + ## Implementation + The result of the unfolding is calculated as follows: + + Lsquared = L#L regularisation conditions squared + + epsilon_j = sum_i A_ij vector of efficiencies + + E^-1 = ((A# Vyy^-1 A)+tau^2 Lsquared) + + x = E (A# Vyy^-1 y + tau^2 Lsquared x0 +lambda/2 * epsilon) is the result + + The derivatives + dx_k/dy_i + dx_k/dA_ij + dx_k/d(tau^2) + are calculated for further usage. + + The covariance matrix V_xx is calculated as: + Vxx_ij = sum_kl dx_i/dy_k Vyy_kl dx_j/dy_l + + ## Warning: + The algorithm is based on "standard" matrix inversion, with the + known limitations in numerical accuracy and computing cost for + matrices with large dimensions. + + Thus the algorithm should not used for large dimensions of x and y + nx should not be much larger than 200 + ny should not be much larger than 1000 + +## Proper choice of tau + One of the difficult questions is about the choice of tau. + The method implemented in TUnfold is the L-curve method: + a two-dimensional curve is plotted + x-axis: log10(chisquare) + y-axis: log10(regularisation condition) + In many cases this curve has an L-shape. The best choice of tau is in the + kink of the L + + Within TUnfold a simple version of the L-curve analysis is available. + It tests a given number of points in a predefined tau-range and searches + for the maximum of the curvature in the L-curve (kink position). + if no tau range is given, the range of the scan is determined automatically + + A nice overview of the L-curve method is given in: + The L-curve and Its Use in the Numerical Treatment of Inverse Problems + (2000) by P. C. Hansen, in Computational Inverse Problems in + Electrocardiology, ed. P. Johnston, + Advances in Computational Bioengineering + http://www.imm.dtu.dk/~pch/TR/Lcurve.ps + + ## Alternative Regularisation conditions + Regularisation is needed for most unfolding problems, in order to avoid + large oscillations and large correlations on the output bins. + It means that some extra conditions are applied on the output bins + + Within TUnfold these conditions are posed on the difference (x-x0), where + x: unfolding output + x0: the bias distribution, by default calculated from + the input matrix A. There is a method SetBias() to change the + bias distribution. + The 3rd argument to DoUnfold() is a scale factor applied to the bias + bias_default[j] = sum_i A[i][j] + x0[j] = scaleBias*bias[j] + The scale factor can be used to + (a) completely suppress the bias by setting it to zero + (b) compensate differences in the normalisation between data + and Monte Carlo + + If the regularisation is strong, i.e. large parameter tau, + then the distribution x or its derivatives will look like the bias + distribution. If the parameter tau is small, the distribution x is + independent of the bias. + + Three basic types of regularisation are implemented in TUnfold + + condition | regularisation + -----------------|------------------------------------ + kRegModeNone | none + kRegModeSize | minimize the size of (x-x0) + kRegModeDerivative | minimize the 1st derivative of (x-x0) + kRegModeCurvature | minimize the 2nd derivative of (x-x0) + + kRegModeSize is the regularisation scheme which often is found in + literature. The second derivative is often named curvature. + Sometimes the bias is not discussed, equivalent to a bias scale factor + of zero. + + The regularisation schemes kRegModeDerivative and + kRegModeCurvature have the nice feature that they create correlations + between x-bins, whereas the non-regularized unfolding tends to create + negative correlations between bins. For these regularisation schemes the + parameter tau could be tuned such that the correlations are smallest, + as an alternative to the L-curve method. + + If kRegModeSize is chosen or if x is a smooth function through all bins, + the regularisation condition can be set on all bins together by giving + the appropriate argument in the constructor (see examples above). + + If x is composed of independent groups of bins (for example, + signal and background binning in two variables), it may be necessary to + set regularisation conditions for the individual groups of bins. + In this case, give kRegModeNone in the constructor and specify + the bin grouping with calls to + RegularizeBins() specify a 1-dimensional group of bins + RegularizeBins2D() specify a 2-dimensional group of bins + + Note, the class TUnfoldDensity provides an automatic setup of complex + regularisation schemes + + For ultimate flexibility, the regularisation condition can be set on each + bin individually + -> give kRegModeNone in the constructor and use + RegularizeSize() regularize one bin + RegularizeDerivative() regularize the slope given by two bins + RegularizeCurvature() regularize the curvature given by three bins + AddRegularisationCondition() + define an arbitrary regulatisation condition +*/ /* This file is part of TUnfold. diff --git a/hist/hist/src/TUnfoldBinning.cxx b/hist/hist/src/TUnfoldBinning.cxx index 86ed4191ce4..120337f0446 100644 --- a/hist/hist/src/TUnfoldBinning.cxx +++ b/hist/hist/src/TUnfoldBinning.cxx @@ -6,37 +6,34 @@ // History: // Version 17.0, initial version, numbered in parallel to TUnfold -////////////////////////////////////////////////////////////////////////// -// -// TUnfoldBinning -// -// This class serves as a container of analysis bins -// analysis bins are specified by defining the axes of a distribution. -// It is also possible to have unconnected analysis bins without axis. -// Multiple TUnfoldBinning objects may be arranged in a tree, -// such that a full tree structure of histograms and bins is supported -// -// If you use this software, please consider the following citation -// S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201] -// -// More documentation and updates are available on -// http://www.desy.de/~sschmitt -// -// Functionality -// -// The class gives access to all analysis bins numbered in sequence. -// Such a sequence of bins may be stored in a 1-dimension histogram. -// Correlations between two TUnfoldBinning objects may be stored in -// a 2-dimensional histogram. This type of ordering is required for -// the TUnfold class. -// -// In addition, it is possible to have root histograms, using the -// axes as defined with the distributions. Underflow/overflow bins -// can be included or excluded when mapping bins on root histograms. -// In addition, it is possible to collapse one of more axes when going -// from a N-dimensional distribution to a root histogram. -// -////////////////////////////////////////////////////////////////////////// +/** \class TUnfoldBinning + \ingroup Hist + This class serves as a container of analysis bins + analysis bins are specified by defining the axes of a distribution. + It is also possible to have unconnected analysis bins without axis. + Multiple TUnfoldBinning objects may be arranged in a tree, + such that a full tree structure of histograms and bins is supported + + If you use this software, please consider the following citation + S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201] + + More documentation and updates are available on + http://www.desy.de/~sschmitt + + Functionality + + The class gives access to all analysis bins numbered in sequence. + Such a sequence of bins may be stored in a 1-dimension histogram. + Correlations between two TUnfoldBinning objects may be stored in + a 2-dimensional histogram. This type of ordering is required for + the TUnfold class. + + In addition, it is possible to have root histograms, using the + axes as defined with the distributions. Underflow/overflow bins + can be included or excluded when mapping bins on root histograms. + In addition, it is possible to collapse one of more axes when going + from a N-dimensional distribution to a root histogram. +*/ /* This file is part of TUnfold. diff --git a/hist/hist/src/TUnfoldDensity.cxx b/hist/hist/src/TUnfoldDensity.cxx index 975f90ff423..18c58726c0a 100644 --- a/hist/hist/src/TUnfoldDensity.cxx +++ b/hist/hist/src/TUnfoldDensity.cxx @@ -6,110 +6,104 @@ // History: // Version 17.0, support for density regularisation, complex binning schemes, tau scan -////////////////////////////////////////////////////////////////////////// -// -// TUnfoldDensity : public TUnfoldSys : public TUnfold -// -// TUnfold is used to decompose a measurement y into several sources x -// given the measurement uncertainties and a matrix of migrations A -// -// More details are described with the documentation of TUnfold. -// -// For most applications, it is best to use TUnfoldDensity -// instead of using TUnfoldSys or TUnfold -// -// If you use this software, please consider the following citation -// S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201] -// -// More documentation and updates are available on -// http://www.desy.de/~sschmitt -// -// -// As compared to TUnfold, TUndolfDensity adds the following functionality -// * background subtraction (see documentation of TUnfoldSys) -// * error propagation (see documentation of TUnfoldSys) -// * regularisation schemes respecting the bin widths -// * support for complex, multidimensional input distributions -// -// Complex binning schemes are imposed on the measurements y and -// on the result vector x with the help of the class TUnfoldBinning -// The components of x or y are part of multi-dimensional distributions. -// The bin widths along the relevant directions in these distributions -// are used to calculate bin densities (number of events divided by bin width) -// or to calculate derivatives taking into account the proper distance of -// adjacent bin centers -// -// Complex binning schemes -// ======================= -// in literature on unfolding, the "standard" test case is a -// one-dimensional distribution without underflow or overflow bins. -// The migration matrix is almost diagonal. -// -// This "standard" case is rarely realized for real problems. -// -// Often one has to deal with multi-dimensional input distributions. -// In addition, there are underflow and overflow bins -// or other background bins, possibly determined with the help of auxillary -// measurements -// -// In TUnfoldDensity, such complex binning schemes are handled with the help -// of the class TUnfoldBinning. For each vector there is a tree -// structure. The tree nodes hold multi-dimensiopnal distributions -// -// For example, the "measurement" tree could have two leaves, one for -// the primary distribution and one for auxillary measurements -// -// Similarly, the "truth" tree could have two leaves, one for the -// signal and one for the background. -// -// each of the leaves may then have a multi-dimensional distribution. -// -// The class TUnfoldBinning takes care to map all bins of the -// "measurement" to the one-dimensional vector y. -// Similarly, the "truth" bins are mapped to the vector x. -// -// Choice of the regularisation -// ============================ -// In TUnfoldDensity, two methods are implemented to determine tau**2 -// (1) ScanLcurve() locate the tau where the L-curve plot has a "kink" -// this function is implemented in the TUnfold class -// (2) ScanTau() finds the solution such that some variable -// (e.g. global correlation coefficient) is minimized -// this function is implemented in the TUnfoldDensity class, -// such that the variable could be made depend on the binning scheme -// -// Each of the algorithms has its own advantages and disadvantages -// -// The algorithm (1) does not work if the input data are too similar to the -// MC prediction, that is unfolding with tau=0 gives a least-square sum -// of zero. Typical no-go cases of the L-curve scan are: -// (a) the number of measurements is too small (e.g. ny=nx) -// (b) the input data have no statistical fluctuations -// [identical MC events are used to fill the matrix of migrations -// and the vector y] -// -// The algorithm (2) only works if the variable does have a real minimum -// as a function of tau. -// If global correlations are minimized, the situation is as follows: -// The matrix of migration typically introduces negative correlations. -// The area constraint introduces some positive correlation. -// Regularisation on the "size" introduces no correlation. -// Regularisation on 1st or 2nd derivatives adds positive correlations. -// For this reason, "size" regularisation does not work well with -// the tau-scan: the higher tau, the smaller rho, but there is no minimum. -// In contrast, the tau-scan is expected to work well with 1st or 2nd -// derivative regularisation, because at some point the negative -// correlations from migrations are approximately cancelled by the -// positive correlations from the regularisation conditions. -// -// whichever algorithm is used, the output has to be checked: -// (1) The L-curve should have approximate L-shape -// and the final choice of tau should not be at the very edge of the -// scanned region -// (2) The scan result should have a well-defined minimum and the -// final choice of tau should sit right in the minimum -// -//////////////////////////////////////////////////////////////////////////// +/** \class TUnfoldBinning + \ingroup Hist + TUnfold is used to decompose a measurement y into several sources x + given the measurement uncertainties and a matrix of migrations A + + More details are described with the documentation of TUnfold. + + For most applications, it is best to use TUnfoldDensity + instead of using TUnfoldSys or TUnfold + + If you use this software, please consider the following citation + S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201] + + More documentation and updates are available on + http://www.desy.de/~sschmitt + + As compared to TUnfold, TUndolfDensity adds the following functionality + * background subtraction (see documentation of TUnfoldSys) + * error propagation (see documentation of TUnfoldSys) + * regularisation schemes respecting the bin widths + * support for complex, multidimensional input distributions + + Complex binning schemes are imposed on the measurements y and + on the result vector x with the help of the class TUnfoldBinning + The components of x or y are part of multi-dimensional distributions. + The bin widths along the relevant directions in these distributions + are used to calculate bin densities (number of events divided by bin width) + or to calculate derivatives taking into account the proper distance of + adjacent bin centers + + ## Complex binning schemes + in literature on unfolding, the "standard" test case is a + one-dimensional distribution without underflow or overflow bins. + The migration matrix is almost diagonal. + + This "standard" case is rarely realized for real problems. + + Often one has to deal with multi-dimensional input distributions. + In addition, there are underflow and overflow bins + or other background bins, possibly determined with the help of auxillary + measurements + + In TUnfoldDensity, such complex binning schemes are handled with the help + of the class TUnfoldBinning. For each vector there is a tree + structure. The tree nodes hold multi-dimensiopnal distributions + + For example, the "measurement" tree could have two leaves, one for + the primary distribution and one for auxillary measurements + + Similarly, the "truth" tree could have two leaves, one for the + signal and one for the background. + + each of the leaves may then have a multi-dimensional distribution. + + The class TUnfoldBinning takes care to map all bins of the + "measurement" to the one-dimensional vector y. + Similarly, the "truth" bins are mapped to the vector x. + + ## Choice of the regularisation + In TUnfoldDensity, two methods are implemented to determine tau**2 + 1. ScanLcurve() locate the tau where the L-curve plot has a "kink" + this function is implemented in the TUnfold class + 2. ScanTau() finds the solution such that some variable + (e.g. global correlation coefficient) is minimized + this function is implemented in the TUnfoldDensity class, + such that the variable could be made depend on the binning scheme + + Each of the algorithms has its own advantages and disadvantages + + The algorithm (1) does not work if the input data are too similar to the + MC prediction, that is unfolding with tau=0 gives a least-square sum + of zero. Typical no-go cases of the L-curve scan are: + - (a) the number of measurements is too small (e.g. ny=nx) + - (b) the input data have no statistical fluctuations + [identical MC events are used to fill the matrix of migrations + and the vector y] + + The algorithm (2) only works if the variable does have a real minimum + as a function of tau. + If global correlations are minimized, the situation is as follows: + The matrix of migration typically introduces negative correlations. + The area constraint introduces some positive correlation. + Regularisation on the "size" introduces no correlation. + Regularisation on 1st or 2nd derivatives adds positive correlations. + For this reason, "size" regularisation does not work well with + the tau-scan: the higher tau, the smaller rho, but there is no minimum. + In contrast, the tau-scan is expected to work well with 1st or 2nd + derivative regularisation, because at some point the negative + correlations from migrations are approximately cancelled by the + positive correlations from the regularisation conditions. + + whichever algorithm is used, the output has to be checked: + 1. The L-curve should have approximate L-shape + and the final choice of tau should not be at the very edge of the + scanned region + 2. The scan result should have a well-defined minimum and the + final choice of tau should sit right in the minimum +*/ /* This file is part of TUnfold. diff --git a/hist/hist/src/TUnfoldSys.cxx b/hist/hist/src/TUnfoldSys.cxx index 8b21b0533d2..467254335e7 100644 --- a/hist/hist/src/TUnfoldSys.cxx +++ b/hist/hist/src/TUnfoldSys.cxx @@ -11,107 +11,98 @@ // Version 14, remove some print-out, do not add unused sys.errors // Version 13, support for systematic errors -//////////////////////////////////////////////////////////////////////// -// -// TUnfoldSys : public TUnfold -// -// TUnfold is used to decompose a measurement y into several sources x -// given the measurement uncertainties and a matrix of migrations A -// -// TUnfoldSys adds error propagation of systematic errors to TUnfold -// Also, background sources (with errors) can be subtracted. -// -// ************************************************************* -// * For most applications, it is better to use TUnfoldDensity * -// * instead of using TUnfoldSys or TUnfold * -// ************************************************************* -// -// If you use this software, please consider the following citation -// S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201] -// -// More documentation and updates are available on -// http://www.desy.de/~sschmitt -// -// -// The following sources of systematic error are considered in TUnfoldSys -// -// (a) uncorrelated errors on the input matrix histA, taken as the -// errors provided with the histogram. -// These are typically statistical errors from finite Monte Carlo samples -// -// (b) correlated shifts of the input matrix histA. These shifts are taken -// as one-sigma effects when switchig on a given error soure. -// several such error sources may be defined -// -// (c) a systematic error on the regularisation parameter tau -// -// (d) uncorrelated errors on background sources, taken as the errors -// provided with the background histograms -// -// (e) scale errors on background sources -// -// In addition there is the (statistical) uncertainty of the input vector (i) -// -// Source (a) is providede with the original histogram histA -// TUnfoldSys(histA,...) -// -// Sources (b) are added by calls to -// AddSysError() -// -// The systematic uncertainty on tau (c) is set by -// SetTauError() -// -// Backgound sources causing errors of type (d) and (e) are added by -// SubtractBackground() -// -// -// NOTE: -// ====== -// Systematic errors (a), (b), (c) are propagated to the result -// AFTER unfolding -// -// Background errors (d) and (e) are added to the data errors -// BEFORE unfolding -// -// For this reason: -// errors of type (d) and (e) are INCLUDED in the standard error matrix -// and other methods provided by the base class TUnfold: -// GetOutput() -// GetEmatrix() -// ... -// whereas errors of type (a), (b), (c) are NOT INCLUDED in the methods -// provided by the base class TUnfold. -// -// Accessing error matrices: -// =================================== -// The error sources (b),(c) and (e) propagate to shifts of the result. -// These shifts may be accessed as histograms using the methods -// GetDeltaSysSource() corresponds to (b) -// GetDeltaSysTau() corresponds to (c) -// GetDeltaSysBackgroundScale() corresponds to (e) -// The error sources (a) and (d) originate from many uncorrelated errors, -// which in general are NOT uncorrelated on the result vector. -// Thus, there is no corresponding shift of the output vector, only error -// matrices are available -// -// Method to get error matrix corresponds to error sources -// =============================================================== -// GetEmatrixSysUncorr() (a) -// GetEmatrixSysSource() (b) -// GetEmatrixSysTau() (c) -// GetEmatrixSysBackgroundUncorr() (d) -// GetEmatrixSysBackgroundScale() (e) -// GetEmatrixInput() (i) -// GetEmatrix() (i)+(d)+(e) -// GetEmatrixTotal() (i)+(a)+(b)+(c)+(d)+(e) -// -// Error matrices can be added to existing histograms. -// This is useful to retreive the sum of several error matrices. -// If the last argument of the GetEmatrixXXX() methods is set to kFALSE, -// the histogram is not cleared, but the error matrix is simply added to the -// existing histogram -// -//////////////////////////////////////////////////////////////////////// +/** \class TUnfoldSys + \ingroup Hist + TUnfold is used to decompose a measurement y into several sources x + given the measurement uncertainties and a matrix of migrations A + + TUnfoldSys adds error propagation of systematic errors to TUnfold + Also, background sources (with errors) can be subtracted. + + **For most applications, it is better to use TUnfoldDensity + instead of using TUnfoldSys or TUnfold** + + If you use this software, please consider the following citation + S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201] + + More documentation and updates are available on + http://www.desy.de/~sschmitt + + The following sources of systematic error are considered in TUnfoldSys + + (a) uncorrelated errors on the input matrix histA, taken as the + errors provided with the histogram. + These are typically statistical errors from finite Monte Carlo samples + + (b) correlated shifts of the input matrix histA. These shifts are taken + as one-sigma effects when switchig on a given error soure. + several such error sources may be defined + + (c) a systematic error on the regularisation parameter tau + + (d) uncorrelated errors on background sources, taken as the errors + provided with the background histograms + + (e) scale errors on background sources + + In addition there is the (statistical) uncertainty of the input vector (i) + + Source (a) is providede with the original histogram histA + TUnfoldSys(histA,...) + + Sources (b) are added by calls to + AddSysError() + + The systematic uncertainty on tau (c) is set by + SetTauError() + + Backgound sources causing errors of type (d) and (e) are added by + SubtractBackground() + + NOTE: + Systematic errors (a), (b), (c) are propagated to the result + AFTER unfolding + + Background errors (d) and (e) are added to the data errors + BEFORE unfolding + + For this reason: + errors of type (d) and (e) are INCLUDED in the standard error matrix + and other methods provided by the base class TUnfold: + GetOutput() + GetEmatrix() + ... + whereas errors of type (a), (b), (c) are NOT INCLUDED in the methods + provided by the base class TUnfold. + + ## Accessing error matrices: + The error sources (b),(c) and (e) propagate to shifts of the result. + These shifts may be accessed as histograms using the methods + GetDeltaSysSource() corresponds to (b) + GetDeltaSysTau() corresponds to (c) + GetDeltaSysBackgroundScale() corresponds to (e) + The error sources (a) and (d) originate from many uncorrelated errors, + which in general are NOT uncorrelated on the result vector. + Thus, there is no corresponding shift of the output vector, only error + matrices are available + + Method to get error matrix | corresponds to error sources + -----------------------------|--------------------------------- + GetEmatrixSysUncorr() | (a) + GetEmatrixSysSource() | (b) + GetEmatrixSysTau() | (c) + GetEmatrixSysBackgroundUncorr() | (d) + GetEmatrixSysBackgroundScale() | (e) + GetEmatrixInput() | (i) + GetEmatrix() | (i)+(d)+(e) + GetEmatrixTotal() | (i)+(a)+(b)+(c)+(d)+(e) + + Error matrices can be added to existing histograms. + This is useful to retreive the sum of several error matrices. + If the last argument of the GetEmatrixXXX() methods is set to kFALSE, + the histogram is not cleared, but the error matrix is simply added to the + existing histogram +*/ /* This file is part of TUnfold. diff --git a/hist/hist/src/TVirtualFitter.cxx b/hist/hist/src/TVirtualFitter.cxx index cea2f5e3199..bc679ca1710 100644 --- a/hist/hist/src/TVirtualFitter.cxx +++ b/hist/hist/src/TVirtualFitter.cxx @@ -10,9 +10,10 @@ *************************************************************************/ -//______________________________________________________________________________ -// -// Abstract Base Class for Fitting +/** \class TVirtualFitter + \ingroup Hist + Abstract Base Class for Fitting +*/ #include "TROOT.h" #include "TVirtualFitter.h" diff --git a/hist/hist/src/TVirtualGraphPainter.cxx b/hist/hist/src/TVirtualGraphPainter.cxx index ec15324a8d3..b91615b9ab2 100644 --- a/hist/hist/src/TVirtualGraphPainter.cxx +++ b/hist/hist/src/TVirtualGraphPainter.cxx @@ -18,11 +18,10 @@ TVirtualGraphPainter *TVirtualGraphPainter::fgPainter = 0; ClassImp(TVirtualGraphPainter) -//______________________________________________________________________________ -// -// TVirtualGraphPainter is an abstract interface to a histogram painter. -// - +/** \class TVirtualGraphPainter + \ingroup Hist + Abstract interface to a histogram painter +*/ //////////////////////////////////////////////////////////////////////////////// /// Static function returning a pointer to the current graph painter. diff --git a/hist/hist/src/TVirtualHistPainter.cxx b/hist/hist/src/TVirtualHistPainter.cxx index 6e952e43f9a..f8ce552bc58 100644 --- a/hist/hist/src/TVirtualHistPainter.cxx +++ b/hist/hist/src/TVirtualHistPainter.cxx @@ -18,10 +18,10 @@ TClass *TVirtualHistPainter::fgPainter = 0; ClassImp(TVirtualHistPainter) -//______________________________________________________________________________ -// -// TVirtualHistPainter is an abstract interface to a histogram painter. -// +/** \class TVirtualHistPainter + \ingroup Hist + Abstract interface to a histogram painter +*/ //////////////////////////////////////////////////////////////////////////////// -- GitLab