- Timestamp:
- 02/19/03 13:47:29 (22 years ago)
- Location:
- trunk/MagicSoft/Mars
- Files:
-
- 4 added
- 9 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/manalysis/MHillasCalc.cc
r1542 r1762 84 84 } 85 85 86 fHillas = (MHillas*)pList->FindCreateObj( fHilName);86 fHillas = (MHillas*)pList->FindCreateObj("MHillas",fHilName); 87 87 if (!fHillas) 88 88 return kFALSE; -
trunk/MagicSoft/Mars/manalysis/MHillasSrcCalc.cc
r1542 r1762 56 56 // The default is "MHillasSrc" 57 57 // 58 //MHillasSrcCalc::MHillasSrcCalc(const char *src, const char *hil, 59 // const char *name, const char *title) 60 // : fHillas(NULL), fSrcPos(NULL), fHillasSrc(NULL) 61 //{ 62 // fName = name ? name : gsDefName.Data(); 63 // fTitle = title ? title : gsDefTitle.Data(); 64 65 // fSrcName = src; 66 // fHillasName = hil; 67 // fHillasInput = "MHillas"; 68 //} 69 // ------------------------------------------------------------------------- 70 // 58 71 MHillasSrcCalc::MHillasSrcCalc(const char *src, const char *hil, 59 const char *name, const char *title) 72 const char *name, const char *title, 73 const char *hilinput) 60 74 : fHillas(NULL), fSrcPos(NULL), fHillasSrc(NULL) 61 75 { … … 65 79 fSrcName = src; 66 80 fHillasName = hil; 81 fHillasInput = hilinput; 67 82 } 68 83 … … 71 86 Bool_t MHillasSrcCalc::PreProcess(MParList *pList) 72 87 { 73 fHillas = (MHillas*)pList->FindObject( "MHillas");88 fHillas = (MHillas*)pList->FindObject(fHillasInput); 74 89 if (!fHillas) 75 90 { -
trunk/MagicSoft/Mars/manalysis/MHillasSrcCalc.h
r1540 r1762 19 19 TString fSrcName; 20 20 TString fHillasName; 21 TString fHillasInput; 21 22 22 23 Int_t fErrors; … … 29 30 30 31 public: 32 // MHillasSrcCalc(const char *src="MSrcPosCam", const char *hil="MHillasSrc", 33 // const char *name=NULL, const char *title=NULL); 34 31 35 MHillasSrcCalc(const char *src="MSrcPosCam", const char *hil="MHillasSrc", 32 const char *name=NULL, const char *title=NULL); 36 const char *name=NULL, const char *title=NULL, 37 const char *hilinput="MHillas"); 33 38 34 39 ClassDef(MHillasSrcCalc, 1) // task to calculate the source position depandant hillas parameters -
trunk/MagicSoft/Mars/mfilter/MFEventSelector.cc
r1668 r1762 94 94 // the class description above. 95 95 // 96 MFEventSelector::MFEventSelector(const char *name, const char *title) 96 MFEventSelector::MFEventSelector(const char *name, const char *title, 97 const char *read) 97 98 : fNumTotalEvts(-1), fNumSelectEvts(-1), fSelRatio(-1), fNumSelectedEvts(0) 98 99 { 99 100 fName = name ? name : gsDefName.Data(); 100 101 fTitle = title ? title : gsDefTitle.Data(); 102 103 fRead = read; 101 104 } 102 105 … … 127 130 return kFALSE; 128 131 } 129 MRead *read = (MRead*)tlist->FindObject( "MRead");132 MRead *read = (MRead*)tlist->FindObject(fRead); 130 133 if (!read) 131 134 { … … 134 137 } 135 138 fNumTotalEvts = read->GetEntries(); 139 140 *fLog << "MFEventSelector::PreProcess; fNumTotalEvts = " 141 << fNumTotalEvts << endl; 136 142 137 143 SetBit(kNumTotalFromFile); -
trunk/MagicSoft/Mars/mfilter/MFEventSelector.h
r1668 r1762 24 24 25 25 Bool_t fResult; 26 TString fRead; 26 27 27 28 void StreamPrimitive(ofstream &out) const; … … 38 39 public: 39 40 // MFEventSelector(); 40 MFEventSelector(const char *name=NULL, const char *title=NULL); 41 MFEventSelector(const char *name=NULL, const char *title=NULL, 42 const char *read="MRead"); 41 43 ~MFEventSelector(); 42 44 -
trunk/MagicSoft/Mars/mhist/MHHillasExt.cc
r1715 r1762 57 57 // Setup four histograms for Width, Length 58 58 // 59 MHHillasExt::MHHillasExt(const char *name, const char *title) 59 MHHillasExt::MHHillasExt(const char *name, const char *title, 60 const char *hil) 60 61 : fMm2Deg(1), fUseMmScale(kTRUE) 61 62 { … … 65 66 fName = name ? name : "MHHillasExt"; 66 67 fTitle = title ? title : "Container for extended Hillas histograms"; 68 69 fHilName = hil; 70 //*fLog << "MHHillasExt : fHilName = " << fHilName << endl; 67 71 68 72 // … … 139 143 Bool_t MHHillasExt::SetupFill(const MParList *plist) 140 144 { 141 TObject *obj = plist->FindObject( "MHillas");145 TObject *obj = plist->FindObject(fHilName, "MHillas"); 142 146 if (!obj) 143 147 { 144 *fLog << err << dbginf << "Sorry 'MHillas' not found in parameter list... aborting." << endl; 148 *fLog << err << dbginf << "Sorry '" << fHilName 149 << "' not found in parameter list... aborting." << endl; 145 150 return kFALSE; 146 151 } -
trunk/MagicSoft/Mars/mhist/MHHillasExt.h
r1574 r1762 25 25 Bool_t fUseMmScale; 26 26 27 TString fHilName; 28 27 29 public: 28 MHHillasExt(const char *name=NULL, const char *title=NULL); 30 //MHHillasExt(const char *name=NULL, const char *title=NULL); 31 MHHillasExt(const char *name=NULL, const char *title=NULL, 32 const char *hil="MHillas"); 33 29 34 ~MHHillasExt(); 30 35 -
trunk/MagicSoft/Mars/mhist/MHMatrix.cc
r1711 r1762 16 16 ! 17 17 ! 18 ! Author(s): Thomas Bretz 2002 <mailto:tbretz@astro.uni-wuerzburg.de> 18 ! Author(s): Thomas Bretz 2002 <mailto:tbretz@astro.uni-wuerzburg.de> 19 ! Rudy Boeck 2003 <mailto: 20 ! Wolfgang Wittek2003 <mailto:wittek@mppmu.mpg.de> 19 21 ! 20 ! Copyright: MAGIC Software Development, 2000-200 222 ! Copyright: MAGIC Software Development, 2000-2003 21 23 ! 22 24 ! … … 52 54 #include <TArrayI.h> 53 55 56 #include <TH1.h> 57 #include <TCanvas.h> 58 #include <TRandom3.h> 59 54 60 #include "MLog.h" 55 61 #include "MLogManip.h" … … 134 140 } 135 141 142 // -------------------------------------------------------------------------- 143 // 136 144 void MHMatrix::AddColumns(MDataArray *matrix) 137 145 { … … 341 349 } 342 350 351 // -------------------------------------------------------------------------- 352 // 343 353 const TMatrix *MHMatrix::InvertPosDef() 344 354 { … … 541 551 } 542 552 553 // -------------------------------------------------------------------------- 554 // 543 555 const TArrayI MHMatrix::GetIndexOfSortedColumn(Int_t ncol, Bool_t desc) const 544 556 { … … 558 570 } 559 571 572 // -------------------------------------------------------------------------- 573 // 560 574 void MHMatrix::SortMatrixByColumn(Int_t ncol, Bool_t desc) 561 575 { … … 577 591 } 578 592 593 // -------------------------------------------------------------------------- 594 // 579 595 Bool_t MHMatrix::Fill(MParList *plist, MTask *read) 580 596 { … … 617 633 // -------------------------------------------------------------------------- 618 634 // 619 // Return a comma seperated list of all data members used in the matrix.620 // This is mainly used in MTask::AddToBranchList621 635 // 622 636 void MHMatrix::ReduceNumberOfRows(UInt_t numrows, const TString opt) … … 674 688 } 675 689 690 // ------------------------------------------------------------------------ 691 // 692 // Define the reference matrix 693 // refcolumn number of the column containing the variable, for which a 694 // target distribution may be given; 695 // if refcolumn is negative the target distribution will be set 696 // equal to the real distribution; the events in the reference 697 // matrix will then be simply a random selection of the events 698 // in the original matrix. 699 // thsh histogram containing the target distribution of the variable 700 // nmaxevts maximum number of events in the reference matrix 701 // 702 Bool_t MHMatrix::DefRefMatrix(const Int_t refcolumn, const TH1F *thsh, 703 const Int_t nmaxevts) 704 { 705 if (!fM.IsValid()) 706 { 707 *fLog << "MHMatrix::DefRefMatrix; matrix not initialized" << endl; 708 return kFALSE; 709 } 710 711 const TH1F *fHthsh; // target distribution 712 TH1F *fHthd; // normalization factors 713 TH1F *fHth; // distribution before 714 TH1F *fHthaft; // distribution after renormalization 715 716 Int_t frefcol; // number of the column (count from 0) containing 717 // the variable for which the target distribution is given 718 Int_t fnbins; // histogram parameters for thsh 719 Float_t ffrombin; // 720 Float_t ftobin; // 721 Float_t fdbin; // 722 723 TArrayF fnormfac; // array of normalization factors 724 TMatrix fMnew; // matrix fM having the requested distribution 725 // and the requested number of rows; 726 // this is the matrix to be used in the g/h separation 727 728 729 //=================================================================== 730 //--------------------------------------------------------- 731 // Calculate normalization factors 732 // 733 734 fHthsh = thsh; 735 736 Int_t fNrows = fM.GetNrows(); 737 Int_t fNcols = fM.GetNcols(); 738 Int_t fColLwb = fM.GetColLwb(); 739 Int_t fRowLwb = fM.GetRowLwb(); 740 741 // print characteristics 742 // 743 //*fLog << "MHMatrix::CalcNormFactors; Event matrix (input): " << endl; 744 //*fLog << " fNrows, fNcols, fColLwb, fRowLwb = " << fNrows << ", " 745 // << fNcols << ", " << fColLwb << ", " << fRowLwb << endl; 746 747 748 //--------------------------------------------------------- 749 // 750 // if refcol < 0 : select reference events randomly 751 // i.e. set the normaliztion factotrs equal to 1.0 752 753 if (refcolumn<0) 754 { 755 frefcol = -refcolumn; 756 } 757 else 758 { 759 frefcol = refcolumn; 760 } 761 762 763 //--------------------------------------------------------- 764 765 if (fHthsh == NULL) 766 { 767 *fLog << "MHMatrix::DefRefMatrix; target distribution has to be defined" 768 << endl; 769 return kFALSE; 770 } 771 772 fnbins = fHthsh->GetNbinsX(); 773 ffrombin = fHthsh->GetBinLowEdge(1); 774 ftobin = fHthsh->GetBinLowEdge(fnbins+1); 775 fdbin = (ftobin-ffrombin)/fnbins; 776 //*fLog << "CalcNormFactors; Reference column (count from 0) = " 777 // << frefcol << endl; 778 //*fLog << "CalcNormFactors; Histo params : " << fnbins << ", " << ffrombin 779 // <<", "<< ftobin << endl; 780 781 782 //--------------------------------------------------------- 783 // ===== set up the real histogram 784 // 785 fHth = new TH1F("th", "data at input", fnbins, ffrombin, ftobin); 786 for (Int_t j=1; j<=fNrows; j++) 787 { 788 fHth->Fill(fM(j-1,frefcol)); 789 //*fLog << "j, frefcol, fM(j-1,frefcol) = " << j << ", " << frefcol 790 // << ", " << fM(j-1,frefcol) << endl; 791 } 792 793 // if refcolumn<0 set target distribution equal to the real distribution 794 // in order to obtain normalization factors = 1.0 795 //if (refcolumn<0) 796 //{ 797 // for (Int_t j=1; j<=fnbins; j++) 798 // { 799 // float cont = fHth-> GetBinContent(j); 800 // fHthsh->SetBinContent(j, cont); 801 // } 802 //} 803 804 //--------------------------------------------------------- 805 // ===== obtain correction factors 806 // 807 fHthd = new TH1F ("thd", "correction factors", fnbins, ffrombin, ftobin); 808 809 Double_t cmax = 0.; 810 Float_t a; 811 Float_t b; 812 Float_t c; 813 for (Int_t j=1; j<=fnbins; j++) 814 { 815 a = fHth->GetBinContent(j); 816 817 // if refcolumn < 0 set the correction factors equal to 1 818 if ( refcolumn>=0.0 ) 819 b = fHthsh->GetBinContent(j); 820 else 821 b = a; 822 823 if ( a<0.0 || b<0.0 ) 824 { 825 *fLog << "MHMatrix::DefRefMatrix; renormalization of input matrix is not possible" << endl; 826 *fLog << " a, b = " << a << ", " << b << endl; 827 return kFALSE; 828 } 829 830 if (a == 0.0) 831 c = 0.0; 832 else 833 c = b/a; 834 835 fHthd->SetBinContent(j, c); 836 if (cmax < c) cmax = c; 837 //*fLog << j << ", " << a << ", " << b << ", " << c << endl; 838 } 839 840 //*fLog << "MHMatrix::CalcNormFactors; maximum ratio between histogram and target " 841 // << cmax << endl; 842 843 if (cmax <= 0.0) 844 { 845 *fLog << "MHMatrix::CalcNormFactors; maximum ratio is LE zero" << endl; 846 return kFALSE; 847 } 848 849 fnormfac.Set(fnbins); 850 Float_t minbin = 1.0; 851 for (Int_t j=1; j<=fnbins; j++) 852 { 853 float c = (fHthd->GetBinContent(j)) / cmax; 854 fnormfac[j-1] = c; 855 if (minbin > c && c>0) minbin = c; 856 fHthd->SetBinContent(j, c); 857 } 858 //*fLog << "MHMatrix::CalcNormFactors; maximum weight = 1, minimum non-zero weight = " 859 // << minbin <<endl; 860 861 862 //--------------------------------------------------------- 863 //=================================================================== 864 865 866 // ==== compress matrix fM into Mnew 867 // 868 // get random access 869 TArrayF ranx(fNrows); 870 TRandom3 *fRnd = new TRandom3(0); 871 872 for (Int_t i=0;i<fNrows;i++) ranx[i] = fRnd->Rndm(i); 873 874 TArrayI ind(fNrows); 875 TMath::Sort(fNrows, ranx.GetArray(), ind.GetArray(),kTRUE); 876 877 //*fLog << "MHMatrix::DefRefMatrix; permuting array " <<endl; 878 //for (Int_t i=0; i<10; i++) cout << ind[i] << " "; 879 //*fLog << endl; 880 881 //--------------------------------------------------------- 882 // go through matrix and keep events depending on content of 883 // reference column 884 // 885 //*fLog << "MHMatrix::DefRefMatrix; pass at most " << nmaxevts 886 // << " events, in random order " << endl; 887 888 889 //--------------------------------------------------------- 890 // define new matrix 891 Int_t evtcount = 0; 892 TMatrix Mnew_tmp (fNrows,fNcols); 893 894 TArrayF cumul_weight(fNrows); // keep track for each bin how many events 895 896 for (Int_t ir=0; ir<fNrows; ir++) cumul_weight[ir]=0; 897 898 //-------------------------------- 899 // select events 900 // 901 fHthaft = new TH1F ("thaft","data after selection", fnbins, ffrombin, ftobin); 902 for (Int_t ir=0; ir<fNrows; ir++) 903 { 904 Int_t indref= (fM(ind[ir],frefcol)-ffrombin)/fdbin; 905 cumul_weight[indref] += fnormfac[indref]; 906 if (cumul_weight[indref]>0.5) 907 { 908 cumul_weight[indref] -= 1.; 909 if (evtcount++ < nmaxevts) 910 { 911 fHthaft->Fill(fM(ind[ir],frefcol)); 912 for (Int_t ic=0; ic<fNcols; ic++) 913 Mnew_tmp(evtcount-1,ic) = fM(ind[ir],ic); 914 } 915 else break; 916 } 917 } 918 evtcount--; 919 //*fLog << "MHMatrix::DefRefMatrix; passed " << evtcount 920 // << " events" << endl; 921 922 //-------------------------------- 923 // reduce size 924 925 TMatrix Mnew(evtcount, fNcols); 926 fM.ResizeTo (evtcount, fNcols); 927 fNumRow = evtcount; 928 929 for (Int_t ir=0;ir<evtcount; ir++) 930 for (Int_t ic=0;ic<fNcols;ic++) 931 { 932 Mnew(ir,ic) = Mnew_tmp(ir,ic); 933 fM(ir,ic) = Mnew_tmp(ir,ic); 934 } 935 936 if (evtcount < nmaxevts) 937 { 938 *fLog << "MHMatrix::DefRefMatrix; Warning : The reference sample contains less events (" << evtcount << ") than required (" << nmaxevts << ")" << endl; 939 } 940 941 942 //--------------------------------------------------------- 943 // test: print new matrix (part) 944 *fLog << "MHMatrix::DefRefMatrix; Event matrix (output) :" << endl; 945 *fLog << "MHMatrix::DefRefMatrix; Nrows, Ncols = " << Mnew.GetNrows() 946 << " " << Mnew.GetNcols() << endl; 947 for (Int_t ir=0;ir<10; ir++) 948 { 949 *fLog <<ir <<" "; 950 for (Int_t ic=0;ic<13;ic++) cout<<Mnew(ir,ic)<<" "; 951 *fLog <<endl; 952 } 953 954 // test print new bin contents 955 *fLog << "MHMatrix::DefRefMatrix; new histogram: " << endl; 956 for (Int_t j=1; j<=fnbins; j++) 957 { 958 float a = fHthaft->GetBinContent(j); 959 if (a>0) *fLog << j << " "<< a << " "; 960 } 961 *fLog <<endl; 962 963 //--------------------------------------------------------- 964 // ==== plot four histograms 965 TCanvas *th1 = new TCanvas (fName, fName, 1); 966 th1->Divide(2,2); 967 968 th1->cd(1); 969 ((TH1F*)fHthsh)->Draw(); // target 970 971 th1->cd(2); 972 ((TH1F*)fHth)->Draw(); // real histogram before 973 974 th1->cd(3); 975 ((TH1F*)fHthd) -> Draw(); // correction factors 976 977 th1->cd(4); 978 ((TH1F*)fHthaft) -> Draw(); // histogram after 979 980 //--------------------------------------------------------- 981 // --- write onto output file 982 // 983 //TFile *outfile = new TFile(fileNameout, "RECREATE", ""); 984 //Mnew.Write(fileNameout); 985 //outfile->Close(); 986 987 return kTRUE; 988 } 989 676 990 // -------------------------------------------------------------------------- 677 991 // … … 691 1005 692 1006 1007 1008 1009 1010 -
trunk/MagicSoft/Mars/mhist/MHMatrix.h
r1711 r1762 13 13 #endif 14 14 15 #include <TArrayF.h> 16 15 17 class TArrayI; 18 class TArrayF; 16 19 20 class TH1F; 17 21 class MTask; 18 22 class MParList; … … 32 36 33 37 MDataArray *fData; // List of data members (columns) 38 34 39 35 40 enum { … … 94 99 Int_t Read(const char *name); 95 100 101 Bool_t DefRefMatrix (const Int_t refcolumn=-1, const TH1F *thsh=NULL, 102 const Int_t nmaxevts=0); 103 96 104 ClassDef(MHMatrix, 1) // Multidimensional Matrix to store events 97 105 }; … … 99 107 #endif 100 108 109
Note:
See TracChangeset
for help on using the changeset viewer.