Ignore:
Timestamp:
02/03/05 12:42:21 (20 years ago)
Author:
mazin
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mtemp/mmpi
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/mtemp/mmpi/SupercutsONOFFClasses/MFindSupercutsONOFF.cc

    r4411 r6232  
    645645    // The default is not cut, i.e. all values (0-1) are taken
    646646   
    647     fThetaMin = 0; // in miliradians  // FIXME: change name
    648     fThetaMax = 90; // new: in deg; old: (1570) in miliradians // FIXME: change name
    649    
     647    //fThetaMin = 0; // in miliradians  // FIXME: change name
     648    //fThetaMax = 1570; // in miliradians // FIXME: change name
     649   
     650    fThetaMin = 0; // degrees (new in magic)
     651    fThetaMax = 90; // degrees (new in magic)
     652
     653
    650654
    651655    fAlphaSig = 20; // By default, signal is expected in alpha<20 for CT1
     
    770774
    771775    // Angular cuts not yet implemented ...
    772     Double_t ConvMMToDeg = 0.00337034;
    773    
    774     fDistCutLow = 0.4/ConvMMToDeg;
    775     fDistCutUp = 1.5/ConvMMToDeg;
    776    
    777     fLengthCutLow = 0.1/ConvMMToDeg;
    778     fLengthCutUp = 1/ConvMMToDeg;
    779 
    780     fWidthCutLow = 0.07/ConvMMToDeg;
    781     fWidthCutUp = 1/ConvMMToDeg;
     776    fConvMMToDeg = 0.00337034;
     777
     778    fDistCutLow = 0.4/fConvMMToDeg;
     779    fDistCutUp = 1.5/fConvMMToDeg;
     780   
     781    fLengthCutLow = 0.1/fConvMMToDeg;
     782    fLengthCutUp = 1/fConvMMToDeg;
     783
     784    fWidthCutLow = 0.07/fConvMMToDeg;
     785    fWidthCutUp = 1/fConvMMToDeg;
    782786
    783787    // ENDTMP
     
    11411145// in previous functions...
    11421146
    1143     fThetaMin = int(ThetaMin*1000.0);
    1144     fThetaMax = int(ThetaMax*1000.0);
     1147  //    fThetaMin = int(ThetaMin*1000.0);
     1148  //  fThetaMax = int(ThetaMax*1000.0);
     1149
     1150
     1151  // in new magic theta is given in deg
     1152  // 0.5 is added so that the rounding to integer is correct
     1153  fThetaMin = int(ThetaMin + 0.5);
     1154  fThetaMax = int(ThetaMax + 0.5);
    11451155
    11461156    fThetaRangeString = ("ThetaRange");
     
    11481158    fThetaRangeString += ("_");
    11491159    fThetaRangeString += (fThetaMax);
    1150     fThetaRangeString += ("mRad");
     1160    // fThetaRangeString += ("mRad");
     1161    fThetaRangeString += ("_Degrees"); // new in magic
     1162
    11511163
    11521164    return kTRUE;
     
    11681180
    11691181
     1182
     1183    return kTRUE;
     1184}
     1185
     1186
     1187Bool_t MFindSupercutsONOFF::SetFilters(Double_t LeakageMax, Double_t DistMax, Double_t DistMin)
     1188{
     1189
     1190    fDistCutLow = DistMin/fConvMMToDeg;
     1191    fDistCutUp  = DistMax/fConvMMToDeg;
     1192    fLeakageMax = LeakageMax;
     1193   
     1194    *fLog << "MFindSupercutsONOFF::SetSizeRange" << endl
     1195          << "Data matrices will be filled with events whose MHillasSrc.fDist " << endl
     1196          << "is in the range "
     1197          << fDistCutLow <<"-"<< fDistCutUp << " degrees" << endl
     1198          << "and fLeakage " << "< " << fLeakageMax << endl;
    11701199
    11711200    return kTRUE;
     
    16821711
    16831712
    1684     // TMP
     1713   
    16851714    // Cuts in Size,
    16861715
     
    16961725    SizeCutMax.SetInverted();
    16971726   
    1698 
    1699 
    1700 
    1701     // ENDTMP
     1727    // Cuts in Dist
     1728    TString DistCutMinString ("MHillasSrc.fDist>");
     1729    DistCutMinString += fDistCutLow;
     1730    MContinue DistCutMin(DistCutMinString);
     1731    DistCutMin.SetInverted();
     1732
     1733    TString DistCutMaxString ("MHillasSrc.fDist<");
     1734    DistCutMaxString += fDistCutUp;
     1735    MContinue DistCutMax(DistCutMaxString);
     1736    DistCutMax.SetInverted();
     1737
     1738
     1739    // Cuts in leakage
     1740    TString LeakCutMaxString ("MNewImagePar.fLeakage1<");
     1741    LeakCutMaxString += fLeakageMax;
     1742    MContinue LeakCutMax(LeakCutMaxString);
     1743    LeakCutMax.SetInverted();
    17021744
    17031745
     
    17061748    Double_t prob = whichfractiontrain /(whichfractiontrain+whichfractiontest);
    17071749                   
     1750   
     1751   
     1752
    17081753    MFRandomSplit split(prob);
    17091754
     
    17351780
    17361781    tlist.AddToList(&read);
    1737 //    tlist.AddToList(&ThetaCutMin);
    1738 //    tlist.AddToList(&ThetaCutMax);
    1739    
    1740     //TMP
     1782    tlist.AddToList(&ThetaCutMin);
     1783    tlist.AddToList(&ThetaCutMax);
     1784   
     1785   
    17411786    tlist.AddToList(&SizeCutMin);
    17421787    tlist.AddToList(&SizeCutMax);
    17431788
    1744     //ENDTMP
     1789    tlist.AddToList(&DistCutMin);
     1790    tlist.AddToList(&DistCutMax);
     1791    tlist.AddToList(&LeakCutMax);
     1792
     1793   
    17451794    tlist.AddToList(&split);
    17461795    tlist.AddToList(&filltrain);
     
    18481897
    18491898
    1850     // TMP
     1899   
    18511900    // Cuts in Size,
    18521901
     
    18631912   
    18641913
    1865 
    1866 
    1867     // ENDTMP
     1914    // Cuts in Dist
     1915    TString DistCutMinString ("MHillasSrc.fDist>");
     1916    DistCutMinString += fDistCutLow;
     1917    MContinue DistCutMin(DistCutMinString);
     1918    DistCutMin.SetInverted();
     1919
     1920    TString DistCutMaxString ("MHillasSrc.fDist<");
     1921    DistCutMaxString += fDistCutUp;
     1922    MContinue DistCutMax(DistCutMaxString);
     1923    DistCutMax.SetInverted();
     1924
     1925
     1926    // Cuts in leakage
     1927    TString LeakCutMaxString ("MNewImagePar.fLeakage1<");
     1928    LeakCutMaxString += fLeakageMax;
     1929    MContinue LeakCutMax(LeakCutMaxString);
     1930    LeakCutMax.SetInverted();
     1931
     1932
     1933
     1934
     1935
     1936   
    18681937
    18691938
     
    18711940   
    18721941    Double_t prob = whichfractiontrain/(whichfractiontrain + whichfractiontest);
    1873                    
     1942
     1943   
     1944
    18741945    MFRandomSplit split(prob);
    18751946
     
    19031974
    19041975    tlist.AddToList(&read);
    1905 //    tlist.AddToList(&ThetaCutMin);
    1906 //    tlist.AddToList(&ThetaCutMax);
     1976    tlist.AddToList(&ThetaCutMin);
     1977    tlist.AddToList(&ThetaCutMax);
    19071978   
    19081979       
    1909     //TMP
     1980   
    19101981    tlist.AddToList(&SizeCutMin);
    19111982    tlist.AddToList(&SizeCutMax);
    19121983
    1913     //ENDTMP
     1984    tlist.AddToList(&DistCutMin);
     1985    tlist.AddToList(&DistCutMax);
     1986    tlist.AddToList(&LeakCutMax);
     1987
     1988   
    19141989
    19151990    tlist.AddToList(&split);
     
    19752050
    19762051
     2052
     2053/// **********///
     2054
     2055// --------------------------------------------------------------------------
     2056//
     2057// Read only training matrices ON and OFF
     2058//
     2059//
     2060
     2061Bool_t MFindSupercutsONOFF::ReadMatrixTrain(const TString &filetrainON, const TString &filetrainOFF)
     2062{
     2063  //--------------------------
     2064  // read in training matrix ON sample
     2065
     2066  TFile filetr(filetrainON);
     2067  fMatrixTrain->Read("MatrixTrain");
     2068  fMatrixTrain->Print("SizeCols");
     2069
     2070  *fLog << "MFindSupercuts::ReadMatrixTrain; Training matrix for ON sample was read in from file '"
     2071        << filetrainON << "'" << endl;
     2072  filetr.Close();
     2073
     2074
     2075  // read in training matrix OFF sample
     2076
     2077  TFile filetrOFF(filetrainOFF);
     2078  fMatrixTrainOFF->Read("MatrixTrainOFF");
     2079  fMatrixTrainOFF->Print("SizeCols");
     2080
     2081  *fLog << "MFindSupercutsONOFF::ReadMatrixTrain; Training matrix for OFF sample was read in from file '"
     2082        << filetrainOFF << "'" << endl;
     2083  filetrOFF.Close();
     2084
     2085  return kTRUE; 
     2086
     2087}
     2088
     2089// --------------------------------------------------------------------------
     2090//
     2091// Read only test matrices ON and OFF
     2092//
     2093//
     2094
     2095Bool_t MFindSupercutsONOFF::ReadMatrixTest(const TString &filetestON, const TString &filetestOFF)
     2096{
     2097  //--------------------------
     2098  // read in testing matrix ON sample
     2099
     2100//--------------------------
     2101  // read in test matrix for ON sample
     2102 
     2103  TFile filete(filetestON);
     2104  fMatrixTest->Read("MatrixTest");
     2105  fMatrixTest->Print("SizeCols");
     2106
     2107  *fLog << "MFindSupercuts::ReadMatrixTest; Test matrix for ON sample was read in from file '"
     2108        << filetestON << "'" << endl;
     2109  filete.Close();
     2110 
     2111
     2112   //--------------------------
     2113  // read in test matrix for OFF sample
     2114 
     2115  TFile fileteOFF(filetestOFF);
     2116  fMatrixTestOFF->Read("MatrixTestOFF");
     2117  fMatrixTestOFF->Print("SizeCols");
     2118
     2119  *fLog << "MFindSupercutsONOFF::ReadMatrixTest; Test matrix for OFF sample was read in from file '"
     2120        << filetestOFF << "'" << endl;
     2121  filete.Close();
     2122
     2123
     2124
     2125  return kTRUE; 
     2126
     2127}
     2128
     2129
    19772130/// **********///
    19782131
     
    19992152  //--------------------------
    20002153  // read in test matrix
    2001 
     2154 
    20022155  TFile filete(filetest);
    20032156  fMatrixTest->Read("MatrixTest");
     
    20072160        << filetest << "'" << endl;
    20082161  filete.Close();
    2009 
     2162 
    20102163  return kTRUE; 
    20112164}
     
    20322185  //--------------------------
    20332186  // read in test matrix
    2034 
     2187 
    20352188  TFile filete(filetest);
    20362189  fMatrixTestOFF->Read("MatrixTestOFF");
     
    20402193        << filetest << "'" << endl;
    20412194  filete.Close();
    2042 
     2195 
    20432196  return kTRUE; 
    20442197}
  • trunk/MagicSoft/Mars/mtemp/mmpi/SupercutsONOFFClasses/MFindSupercutsONOFF.h

    r4411 r6232  
    5959    Double_t  fWidthCutUp;
    6060
     61    Double_t  fLeakageMax;
     62
     63    Double_t fConvMMToDeg;
    6164
    6265    // ENDTMP
     
    8588  Double_t fAlphaBkgMax;
    8689 
    87   Int_t fThetaMin; // Cuts in ThetaOrig.fVal (in mili radians!!!)
    88   Int_t fThetaMax; // Cuts in ThetaOrig.fVal (in mili radians !!!)
     90  Int_t fThetaMin; // Cut in Theta value (in degrees)
     91  Int_t fThetaMax; // Cuts in Theta value (in degrees)
    8992  TString fThetaRangeString;
    9093
     
    314317  Bool_t CheckAlphaSigBkg();
    315318
     319
     320
     321
     322
    316323  void SetUseFittedQuantities (Bool_t b)
    317324      {fUseFittedQuantities = b;}
     
    371378  Bool_t ReadMatrixOFF( const TString &filetrainOFF, const TString &filetestOFF);
    372379
     380
     381  Bool_t ReadMatrixTrain(const TString &filetrainON, const TString &filetrainOFF);
     382  Bool_t ReadMatrixTest(const TString &filetestON, const TString &filetestOFF);
     383
     384
     385
    373386  Bool_t ComputeNormFactorTrain ();
    374387  Bool_t ComputeNormFactorTest ();
     
    395408
    396409  Bool_t SetSizeRange(Double_t SizeMin, Double_t SizeMax);
     410
     411  Bool_t SetFilters(Double_t LeakageMax, Double_t DistMax, Double_t DistMin);
     412
    397413
    398414
  • trunk/MagicSoft/Mars/mtemp/mmpi/SupercutsONOFFClasses/MFindSupercutsONOFFThetaLoop.cc

    r4411 r6232  
    110110    //---------------------------
    111111
    112 
     112    // *********************************
     113    // ** NOT TRUE ANY LONGER... ** //
    113114    // Cuts (low and up) in variable ThetaOrig.fVal
    114115    // The default is not cut, i.e. all values (0-1) are taken
    115116    // fThetaOrig.fVal is measured in Radians; thus 1 = 57 degrees.
    116117
     118    //    fThetaMin = 0.0;
     119    //    fThetaMax = 1.0;
     120    // *********************************
     121
     122
     123    // For MAGIC new (July 2004) Theta is measured in degrees
     124
    117125    fThetaMin = 0.0;
    118     fThetaMax = 1.0;
     126    fThetaMax = 90.0;
     127
     128
    119129
    120130    fAlphaSig = 20; // By default, signal is expected in alpha<20 for CT1
     
    545555
    546556
    547     // Give a "soft warning" if Vector has increasing values
    548     // of CosThetas; which means decreasing values of Thetas
    549     /*
    550     if (fCosThetaRangeVector[0] <= fCosThetaRangeVector[1])
    551     {
    552         *fLog << "MFindSupercutsONOFFThetaLoop::SetCosThetaRangeVector; "
    553               << endl
    554               << "Values in vector fCosThetaRangeVector are in increasing "
    555               << "order, i.e., values in theta will be in decreasing order"
    556               << endl
    557               << "DO NOT forget it !!" << endl;
    558     }
    559 
    560     */
     557   
    561558
    562559    // fCosThetaBinCenterVector is defined
     
    622619    // Now costheta values have to be converted to Radians
    623620    // Function TMath::ACos(x) gives angle in Radians
    624 
     621   
    625622   
    626623    fThetaMin = TMath::ACos(fThetaMin);
    627624    fThetaMax = TMath::ACos(fThetaMax);
     625
     626   
     627    // For new MAGIC (July 2004), angles are given in degrees.
     628    // So I convert angles into degrees by multiplying by a the
     629    // conversion rad/degree
     630
     631    fThetaMin = fThetaMin * TMath::RadToDeg();
     632    fThetaMax = fThetaMax * TMath::RadToDeg();
     633
     634
    628635
    629636    // If fThetaMin > fThetaMax means that fCosThetaRangeVector
     
    644651          << "MFindSupercutsONOFFThetaLoop::SetThetaRange; " << endl
    645652          << "Theta range for this iteration is set to "
    646           << fThetaMin << "-" << fThetaMax << " Radians." << endl
     653          << fThetaMin << "-" << fThetaMax << " degrees." << endl
    647654          << "-------------------------------------------------"<< endl;
    648655   
     
    932939   
    933940
    934     // Histograms references are removed from the current directory
    935 
    936     /* IT DOES NOT WORK !!! THESE STATEMENTS PRODUCE SEGMENTATION FAULT.
    937     fNormFactorTrainHist -> SetDirectory(NULL);
    938    
    939     fNormFactorTestHist -> SetDirectory(NULL);
    940    
    941     fSigmaLiMaTrainHist -> SetDirectory(NULL);
    942    
    943     fSigmaLiMaTestHist -> SetDirectory(NULL);
    944    
    945     fSigmaLiMaTestHist -> SetDirectory(NULL);
    946    
    947     fNexTestHist -> SetDirectory(NULL);
    948     */
     941   
    949942
    950943
     
    10591052    Double_t fofftest)
    10601053{
     1054
     1055  // Check that fractions set by user are within the expected limits
     1056 
     1057  if (fontrain < 0.00 || fontrain > 1.00 ||
     1058      fontest < 0.00 || fontest > 1.00 ||
     1059      fofftrain < 0.00 || fofftrain > 1.00 ||
     1060      fofftest < 0.00 || fontest > 1.00)
     1061    {
     1062      *fLog << err
     1063            << "MFindSupercutsONOFFThetaLoop::SetFractionTrainTestOnOffEvents: "
     1064            << endl
     1065            << "Fractions of ON/OFF data for TRAIN/TEST sample which were specified by user " << endl
     1066            << "are outside the natural limits 0.0-1.0. " << endl
     1067            << "Set reasonable fraction numbers for TRAIN/TEST please. The program will be aborted."
     1068            << endl;
     1069     
     1070      exit(1);
     1071
     1072
     1073    }
     1074 
     1075
     1076  if((fontrain+fontest) > 1.00 ||
     1077     (fofftrain+fofftest) >1.00)
     1078   
     1079    {
     1080      *fLog << err
     1081            << "MFindSupercutsONOFFThetaLoop::SetFractionTrainTestOnOffEvents: "
     1082            << endl
     1083            << "The sum of the fractions of TRAIN/TEST for ON or OFF data which were specified by user " << endl
     1084            << "are larger than 1.0. " << endl
     1085            << "Set reasonable fraction numbers for TRAIN/TEST please. The program will be aborted."
     1086            << endl;
     1087     
     1088      exit(1);
     1089
     1090
     1091    }
     1092
    10611093
    10621094    fWhichFractionTrain = fontrain;
     
    11541186
    11551187    // Path + "OptSCParametersONOFF" + "ThetaRange"
    1156     // + int {(fThetaMin_fThetaMax)*1000} + ".root"
     1188    // + int {(fThetaMin_fThetaMax)} + ".root"
    11571189
    11581190
     
    11601192
    11611193    // Path + "ON(OFF)Data" + "Train(Test)" + "Matrix" + Fraction
    1162     // + "ThetaRange" + int {(fThetaMin_fThetaMax)*1000} + ".root"
     1194    // + "ThetaRange" + int {(fThetaMin_fThetaMax)} + ".root"
    11631195
    11641196
     
    11761208
    11771209
    1178     Int_t ThetaMinTimes1000 = 0;
    1179     Int_t ThetaMaxTimes1000 = 0;
     1210    //    Int_t ThetaMinTimes1000 = 0;
     1211    // Int_t ThetaMaxTimes1000 = 0;
     1212
     1213    // For MAGIC new (July 2004) Theta is given in deg
     1214
     1215    Int_t ThetaMinInt = 0;
     1216    Int_t ThetaMaxInt = 0;
     1217
     1218
    11801219    Int_t FractionONTrainTimes100 = 0;
    11811220    Int_t FractionONTestTimes100 = 0;
     
    12001239       
    12011240       
    1202         ThetaMinTimes1000 = int(fThetaMin*1000);
    1203         ThetaMaxTimes1000 = int(fThetaMax*1000);
    1204         sprintf(ThetaRangeString, "ThetaRange%d_%dmRad",
    1205                 ThetaMinTimes1000, ThetaMaxTimes1000);
     1241        //ThetaMinTimes1000 = int(fThetaMin*1000);
     1242        //ThetaMaxTimes1000 = int(fThetaMax*1000);
     1243        //sprintf(ThetaRangeString, "ThetaRange%d_%dmRad",
     1244        //      ThetaMinTimes1000, ThetaMaxTimes1000);
     1245
     1246        // For New MAGIC (July 2004)
     1247        // Theta is given in deg, no need to multiply by 1000
     1248       
     1249        // 0.5 is added so that the rounding to integer is correct
     1250        ThetaMinInt = int(fThetaMin + 0.5);
     1251        ThetaMaxInt = int(fThetaMax + 0.5);
     1252        sprintf(ThetaRangeString, "ThetaRange%d_%d_Degrees",
     1253                ThetaMinInt, ThetaMaxInt);
     1254
     1255       
    12061256
    12071257        fThetaRangeStringVector[i] = (ThetaRangeString);
     
    12111261        // endtmp
    12121262
    1213         Double_t tmpdouble = 0.0;
    1214         tmpdouble = fWhichFractionTrain*100 - int(fWhichFractionTrain*100);
    1215        
    1216         FractionONTrainTimes100 = tmpdouble < 0.5 ? int(fWhichFractionTrain*100):
    1217             int(fWhichFractionTrain*100) + 1;
    1218        
    1219 
    1220         tmpdouble = fWhichFractionTest*100 - int(fWhichFractionTest*100);
    1221 
    1222         FractionONTestTimes100 = tmpdouble < 0.5 ? int(fWhichFractionTest*100):
    1223             int(fWhichFractionTest*100) + 1;
    1224        
    1225         tmpdouble = fWhichFractionTrainOFF*100 - int(fWhichFractionTrainOFF*100);
    1226 
    1227         FractionOFFTrainTimes100 = tmpdouble < 0.5 ? int(fWhichFractionTrainOFF*100):
    1228             int(fWhichFractionTrainOFF*100) + 1;
    1229 
    1230         tmpdouble = fWhichFractionTestOFF*100 - int(fWhichFractionTestOFF*100);
    1231         FractionOFFTestTimes100 = tmpdouble < 0.5 ? int(fWhichFractionTestOFF*100):
    1232             int(fWhichFractionTestOFF*100) + 1;
     1263        // 0.5 is added so that the rounding to integer is correct
     1264
     1265        FractionONTrainTimes100 = int(fWhichFractionTrain*100 + 0.5);
     1266        FractionONTestTimes100 = int(fWhichFractionTest*100 + 0.5);
     1267        FractionOFFTrainTimes100 = int(fWhichFractionTrainOFF*100 + 0.5);
     1268        FractionOFFTestTimes100 = int(fWhichFractionTestOFF*100 + 0.5);
     1269         
     1270
    12331271       
    12341272       
     
    15251563        FindSupercuts.SetSizeRange(fSizeCutLow, fSizeCutUp);
    15261564
     1565        FindSupercuts.SetFilters(fLeakMax, fDistMax, fDistMin);
    15271566
    15281567
     
    16271666                  << "Reading Data matrices from root files... " << endl;
    16281667
     1668            /*
    16291669            FindSupercuts.ReadMatrix(fTrainMatrixONFilenameVector[i],
    16301670                                     fTestMatrixONFilenameVector[i]);
    16311671            FindSupercuts.ReadMatrixOFF(fTrainMatrixOFFFilenameVector[i],
    16321672                                        fTestMatrixOFFFilenameVector[i]);
     1673            */
     1674
     1675            FindSupercuts.ReadMatrixTrain(fTrainMatrixONFilenameVector[i],
     1676                                          fTrainMatrixOFFFilenameVector[i]);
     1677
     1678            if (fTestParameters)
     1679              {
     1680                FindSupercuts.ReadMatrixTest(fTestMatrixONFilenameVector[i],
     1681                                             fTestMatrixOFFFilenameVector[i]);
     1682              }
    16331683
    16341684        }
     
    23252375       
    23262376            // TEMP   
    2327             // cout << NexVSAlphaNameTmp << endl;
    2328             // cout << SignificanceVSAlphaNameTmp << endl;
     2377            //cout << NexVSAlphaNameTmp << endl;
     2378            //cout << SignificanceVSAlphaNameTmp << endl;
    23292379            // ENDTEMP
    23302380             
     
    23882438
    23892439            AlphaOFFAfterCuts = (TH1F*) rootfile.Get(TmpAlphaOFFHistoName);
     2440
     2441           
     2442
    23902443           
    23912444
     
    28852938{
    28862939
     2940
     2941   *fLog << endl << endl
     2942        << "****************************************************************************" << endl
     2943        << "Combining number of excess events and significance for the  " << endl
     2944        << "selected bins in theta angle using the function " << endl
     2945        << " MFindSupercutsONOFFThetaLoop::ComputeOverallSignificance() " << endl
     2946        << "****************************************************************************" << endl
     2947        << endl;
     2948 
     2949
     2950
     2951 
     2952
     2953
    28872954    // First of all let's check that we have the proper the "drinks" for the party
    28882955
     
    29172984    // Remember that histogram index START AT 1 !!!
    29182985
     2986
     2987
     2988    // Histograms that will contain the alpha distributions
     2989    // for the combined theta bins
     2990
     2991
    29192992    TH1F* CombinedAlphaTrainON = NULL;
    2920     TString CombinedAlphaTrainONName = ("CombinedAlphaTrainON");
     2993    TString CombinedAlphaTrainONName = ("TrainSampleAlphaONCombinedThetaBins");
    29212994    TH1F* CombinedAlphaTrainOFF = NULL;
    2922     TString CombinedAlphaTrainOFFName = ("CombinedAlphaTrainOFF");
     2995    TString CombinedAlphaTrainOFFName = ("TrainSampleAlphaOFFCombinedThetaBins");
    29232996   
    29242997    TH1F* CombinedAlphaTestON = NULL;
    2925     TString CombinedAlphaTestONName = ("CombinedAlphaTestON");
     2998    TString CombinedAlphaTestONName = ("TestSampleAlphaONCombinedThetaBins");
    29262999   
    29273000    TH1F* CombinedAlphaTestOFF = NULL;
    2928     TString CombinedAlphaTestOFFName = ("CombinedAlphaTestOFF");
     3001    TString CombinedAlphaTestOFFName = ("TestSampleAlphaONCombinedThetaBins");
    29293002   
    2930 
    2931    
     3003   
     3004
    29323005
    29333006   
     
    30893162    {
    30903163
     3164
     3165     
     3166      *fLog << endl
     3167        << "****************************************************************************" << endl
     3168        << "Combining data for the TRAIN sample " << endl
     3169        << "****************************************************************************" << endl
     3170        << endl;
     3171
    30913172        // Retrieving NormFactorTrainHisto from root file.
    30923173       
     
    34893570        CombinedAlphaTrainOFF->SetEntries(tmpint);
    34903571
    3491 
     3572        // Check that the number of events in the combined normalized histogram is
     3573        // that of the initial sample
     3574
     3575        /*
    34923576        cout << "SILLY INFO FOR TRAIN SAMPLE: Total number of events Before nomralization and "
    34933577             << " re-re normalized" << endl
    34943578             << TotalNumberOFFTrain << " - " << EventCounter << endl;
    34953579
     3580        */
    34963581    }
    34973582
     
    35063591    if (CombineTestData)
    35073592    {
     3593
     3594      *fLog << endl
     3595        << "****************************************************************************" << endl
     3596        << "Combining data for the TEST sample " << endl
     3597        << "****************************************************************************" << endl
     3598        << endl;
     3599
    35083600
    35093601        // Retrieving NormFactorTestHisto from root file.
     
    39304022
    39314023
     4024        // Check that the number of events in the combined normalized histogram is
     4025        // that of the initial sample
     4026
     4027        /*
    39324028        cout << "SILLY INFO FOR TEST SAMPLE: Total number of events Before nomralization and "
    39334029             << " re-re normalized" << endl
    39344030             << TotalNumberOFFTest << " - " << EventCounter << endl;
     4031
     4032        */
     4033
    39354034
    39364035        /*
     
    39824081    if (SuccessfulThetaBinsVector.GetSize() >= 1)
    39834082    {
     4083
     4084
     4085      *fLog << endl
     4086        << "**************************************************************************************" << endl
     4087        << "Computation of excess events and singnificance for the combined alpha histograms " << endl
     4088        << "**************************************************************************************" << endl
     4089        << endl;
     4090
    39844091        // There is, at least, ONE theta bin for which optimization of
    39854092        // supercuts was SUCCESSFUL, and thus, it is worth to
     
    40264133
    40274134            // Name of psfile where Alpha plot will be stored
    4028             TString psfilename = ("TRAINSampleCombined");
     4135         
     4136            TString psfilename = (fPathForFiles);
     4137            psfilename += ("TRAINSampleCombinedThetaBins");
    40294138           
    40304139            findsig.FindSigmaONOFF(CombinedAlphaTrainON, CombinedAlphaTrainOFF,
     
    40534162
    40544163            // Name of psfile where Alpha plot will be stored
    4055             TString psfilename = ("TESTSampleCombined");
     4164          TString psfilename = (fPathForFiles);
     4165          psfilename += ("TESTSampleCombinedThetaBins");
    40564166           
    40574167            findsig.FindSigmaONOFF(CombinedAlphaTestON, CombinedAlphaTestOFF,
     
    40754185       
    40764186
     4187        *fLog << endl
     4188              << "**************************************************************************************" << endl
     4189              << "Alpha distributions (Train/Test and ON/OFF) for the combined  " << endl
     4190              << " theta bins are stored in the root file " << endl
     4191              << fAlphaDistributionsRootFilename << endl
     4192              << "**************************************************************************************" << endl
     4193        << endl;
     4194
     4195
     4196
     4197        TFile rootfiletowrite(fAlphaDistributionsRootFilename, "UPDATE");
     4198
     4199
     4200       
     4201        // Histograms that will contain the Norm factors for the
     4202        // combined TRAIN/TEST theta bins are created
     4203
     4204       
     4205        TString CombinedNormFactorTrainName = ("TrainSampleNormFactorCombinedThetaBins");
     4206        TH1F CombinedNormFactorTrain (CombinedNormFactorTrainName,
     4207                                      CombinedNormFactorTrainName,
     4208                                      1, 0.0, 1.0);
     4209
     4210        TString CombinedNormFactorTestName = ("TestSampleNormFactorCombinedThetaBins");
     4211        TH1F CombinedNormFactorTest (CombinedNormFactorTestName,
     4212                                     CombinedNormFactorTestName,
     4213                                     1, 0.0, 1.0);
     4214
     4215
     4216        if (CombineTrainData)
     4217        {
     4218
     4219        // Histogram that will contain the Norm factors for the
     4220        // combined TRAIN theta bins are filled
     4221       
     4222       
     4223        CombinedNormFactorTrain.Fill(0.5, MeanNormFactorTrain);
     4224
     4225
     4226        // Train histos are written to root file
     4227       
     4228        CombinedAlphaTrainON -> Write();
     4229        CombinedAlphaTrainOFF -> Write();
     4230        CombinedNormFactorTrain.Write();
     4231
     4232
     4233
     4234        }
     4235
     4236
     4237        if (CombineTestData)
     4238        {
     4239
     4240          // Histogram that will contain the Norm factors for the
     4241        // combined TEST theta bins are filled
     4242
     4243       
     4244          CombinedNormFactorTest.Fill(0.5, MeanNormFactorTest);
     4245
     4246          // Test histos are written to root file
     4247
     4248        CombinedAlphaTestON -> Write();
     4249        CombinedAlphaTestOFF -> Write();
     4250        CombinedNormFactorTest.Write();
     4251       
     4252
     4253        }
     4254
     4255
     4256        // Root file is closed
     4257
     4258        rootfiletowrite.Close();
     4259
    40774260       
    40784261    }
  • trunk/MagicSoft/Mars/mtemp/mmpi/SupercutsONOFFClasses/MFindSupercutsONOFFThetaLoop.h

    r4411 r6232  
    7373  Double_t fSizeCutUp;
    7474
     75  Double_t fLeakMax;
     76  Double_t fDistMax;
     77  Double_t fDistMin;
    7578 
    7679 // Variables for binning of alpha plots
     
    373376   void SetSizeRange(Double_t SizeMin, Double_t SizeMax)
    374377        {fSizeCutLow = SizeMin; fSizeCutUp = SizeMax; }
     378
     379//DM:
     380   void SetFilters(Double_t LeakageMax, Double_t DistMax, Double_t DistMin)
     381        {fLeakMax = LeakageMax; fDistMax = DistMax; fDistMin = DistMin; }
    375382
    376383   // Function that loops over the theta ranges defined by
  • trunk/MagicSoft/Mars/mtemp/mmpi/SupercutsONOFFClasses/MHFindSignificanceONOFF.cc

    r5329 r6232  
    132132
    133133      //-----------------------------
    134       // ignore unwanted points
    135       if (error > 1.e19)
     134      // ignore unwanted points      if (error > 1.e19)
    136135        continue;
    137136
     
    373372
    374373    // allow rebinning of the alpha plot
    375     fRebin = kTRUE;
     374//    fRebin = kTRUE;
     375
     376    fRebin = kFALSE;
    376377
    377378    // allow reducing the degree of the polynomial
     
    396397
    397398    fPrintResultsOntoAlphaPlot = kTRUE;
     399    //fPrintResultsOntoAlphaPlot = kFALSE;
    398400
    399401
     
    872874  if (fDraw)
    873875  {
    874 
    875 
    876876       // TEMPORALLY I will plot fHistOFF and fHistOFFNormalized (with the fits)
    877877
     
    10461046
    10471047      // count bins with zero entry
    1048       if (content <= 0.0)
    1049         {
    1050           fNzeroOFF++;
    1051           // The error of the bin is set to a huge number,
    1052           // so that it does not have any weight in the fit
    1053           fHistOFF->SetBinError(i, dummy);
    1054         }
     1048      if (content < 0.0001)
     1049        fNzeroOFF++;
     1050     
    10551051      // set minimum error
    1056       if (content < 9.0)
     1052      // modified by K. Mase
     1053      if (content < 0.0)
    10571054      {
    10581055        fMlowOFF += 1;
     
    10811078    mean /= ((Double_t) fMbinsOFF);
    10821079    rms  /= ((Double_t) fMbinsOFF);
    1083   }
    1084 
    1085   rms = sqrt( rms - mean*mean );
     1080    rms = sqrt( rms - mean*mean );
     1081  }
     1082
     1083 
    10861084
    10871085  // if there are no events in the background region
     
    11551153
    11561154  fConstantBackg = kFALSE;
    1157 
    1158   // Condition for disabling the fitting procedure and
    1159   // assuming a constant background (before Nov 2004)
    1160 
    1161   // if ( fNzeroOFF > 0  ||  (Double_t)fMlowOFF>0.05*(Double_t)fMbinsOFF )
    1162 
    1163 
    1164   // Condition for disabling the fitting procedure and
    1165   // assuming a constant background (After Nov 01 2004)
    1166   // I softened the condition to allow the fit also in situations
    1167   // where the reduction of the background is such that very
    1168   // few events survived; which is
    1169   // Specially frequent with Random Forest at high Sizes)
    1170 
    1171    if ( (Double_t)fNzeroOFF > 0.1*(Double_t)fMbinsOFF || 
    1172        (Double_t)fMlowOFF > 0.2*(Double_t)fMbinsOFF )
     1155  if ( fNzeroOFF > 0  ||  (Double_t)fMlowOFF>0.05*(Double_t)fMbinsOFF )
    11731156  {
    11741157    *fLog << "MHFindSignificanceONOFF::FitPolynomialOFF; polynomial fit not possible,  fNzeroOFF, fMlowOFF, fMbinsOFF = "
     
    12011184    Double_t val, err;
    12021185    val = mean;
    1203     err = rms; // sqrt( mean / (Double_t)fMbinsOFF );
     1186    err = rms;
     1187   
    12041188
    12051189    fPolyOFF->SetParameter(0, val);
     
    15711555
    15721556      // count bins with zero entry
    1573       if (content <= 0.0)
    1574         {
    1575           fNzero++;
    1576           // The error of the bin is set to a huge number,
    1577           // so that it does not have any weight in the fit
    1578           fHistOFF->SetBinError(i, dummy);
    1579         }
    1580 
     1557      /*if (content <= 0.0)
     1558        fNzero++;*/
     1559     
    15811560      // set minimum error
    1582       if (content < 9.0)
     1561      /*if (content < 9.0)
    15831562      {
    15841563        fMlow += 1;
    15851564        fHist->SetBinError(i, 3.0);
    1586       }
     1565        }*/
    15871566
    15881567      //*fLog << "Take : i, content, error = " << i << ",  "
     
    16811660
    16821661  fConstantBackg = kFALSE;
    1683 
    1684   // Condition for disabling the fitting procedure and
    1685   // assuming a constant background (before Nov 2004)
    1686 
    1687   // if ( fNzero > 0  ||  (Double_t)fMlow>0.05*(Double_t)fMbins )
    1688 
    1689 
    1690   // Condition for disabling the fitting procedure and
    1691   // assuming a constant background (After Nov 01 2004)
    1692   // I softened the condition to allow the fit also in situations
    1693   // where the reduction of the background is such that very
    1694   // few events survived; which is
    1695   // Specially frequent with Random Forest at high Sizes)
    1696 
    1697   if ( (Double_t)fNzero > 0.1*(Double_t)fMbins || 
    1698        (Double_t)fMlow > 0.2*(Double_t)fMbins )
    1699 
     1662  if ( fNzero > 0  ||  (Double_t)fMlow>0.05*(Double_t)fMbins )
    17001663  {
    17011664    *fLog << "MHFindSignificanceONOFF::FitPolynomial; polynomial fit not possible,  fNzero, fMlow, fMbins = "
     
    17281691    Double_t val, err;
    17291692    val = mean;
    1730     err = rms; // sqrt( mean / (Double_t)fMbins );
     1693    err = sqrt( mean / (Double_t)fMbins );
    17311694
    17321695    fPoly->SetParameter(0, val);
     
    27762739     }
    27772740     sum  *= fac*fac;
    2778      
     2741
    27792742     if (sum < 0.0)
    27802743     {
     
    32313194
    32323195    CanvasHistOFFNormalized -> Update();
    3233 
    3234 
    3235     */
     3196*/
    32363197
    32373198    return kTRUE;
     
    32923253            fpoly->DrawCopy("same");
    32933254        }
     3255
     3256        fHistOFFNormalized->SetLineColor(kRed);
     3257        fHistOFFNormalized->DrawCopy("same");
     3258
    32943259    }
    32953260
     
    33163281      {
    33173282        // 10, 1 is white and solid
    3318         fbackg->SetLineColor(10);
     3283        fbackg->SetLineColor(kRed);
    33193284        fbackg->SetLineStyle(1);
    33203285        fbackg->SetLineWidth(4);
    3321         // fbackg->DrawCopy("same"); I do not want to draw it... already too many things.
     3286        //fbackg->DrawCopy("same"); // I do not want to draw it... already too many things.
    33223287      }
    33233288    }
     
    34953460            // fPsFilenameString
    34963461           
    3497             filename += ("AlphaPlotAfterSupercuts.ps");
     3462            filename += ("AlphaPlot.ps");
    34983463            cout << filename << endl;   
    34993464            fCanvas -> SaveAs(filename);
  • trunk/MagicSoft/Mars/mtemp/mmpi/SupercutsONOFFClasses/MSupercutsCalcONOFF.cc

    r4411 r6232  
    133133    fWidthUpperLimit = 0.4;
    134134
     135    // Temp upper limit to get rid of triangular events
     136    fLog10ConcUpperLimit = -0.35;
     137
    135138    fLeakage1UpperLimit = 0.25;
    136139
     
    138141
    139142    fDistLowerLimit = 0.1;
    140     fLengthLowerLimit = 0.09;
    141     fWidthLowerLimit = 0.05;
     143    fLengthLowerLimit = 0.05; // values at October 21th 2004
     144    fWidthLowerLimit = 0.03;  // values at October 21th 2004
    142145
    143146   
     
    237240
    238241}
     242
     243
     244// Function implementing a filter for muon induced images
     245// (Values from Keichi cuts October 18th 2004)
     246
     247// Function returns kTRUE if the event does NOT survive the
     248// muon filter. So kTRUE means event MUST be rejected
     249
     250
     251Bool_t MSupercutsCalcONOFF::MuonFilter (Double_t Size,
     252                                        Double_t Length)
     253{
     254  // Parametrization for cut in length.
     255  // Event is rejected if
     256  // Log10(Length) < k1*(log10(Size)) + k2   (Size<10000photons)
     257  // k1 = 0.286; k2 = -1.91 (October 18th 2004)
     258
     259  // Length < 0.2 deg (Size > 10000)
     260
     261
     262  // Parametrization for cut in Length/Size
     263  // Event is rejected if
     264  // log10(LenOverSize) < k3*(log10(Size)) + k4 (Size<4000photons)
     265  // k3 = -0.780; k4 = -1.71
     266
     267
     268  Double_t LengthUpperLimit = 0.0;
     269  Double_t LengthOverSizeUpperLimit = 0.0;
     270 
     271
     272  Double_t SizeUpperLimitForCutInLength = 10000; // in photons
     273  Double_t SizeUpperLimitForCutInLengthOverSize = 4000; // in photons
     274 
     275  Double_t k1 =  0.286;
     276  Double_t k2 =  -1.91;
     277  Double_t k3 =  -0.78;
     278  Double_t k4 =  -1.71;
     279 
     280
     281 
     282
     283
     284  if (Size <= SizeUpperLimitForCutInLength)
     285    {
     286      LengthUpperLimit = k1 * TMath::Log10(Size) + k2;
     287      LengthUpperLimit = TMath::Power(10, LengthUpperLimit);
     288    }
     289  else
     290    {
     291      LengthUpperLimit = 0.2;
     292    }
     293
     294  // Apply cut in Length
     295 
     296  if (Length <  LengthUpperLimit)
     297        return kTRUE;
     298
     299
     300 
     301  if (Size < SizeUpperLimitForCutInLengthOverSize)
     302    {
     303      LengthOverSizeUpperLimit = k3 * TMath::Log10(Size) + k4;
     304      LengthOverSizeUpperLimit = TMath::Power(10,  LengthOverSizeUpperLimit);
     305 
     306      if (Length/Size <  LengthOverSizeUpperLimit)
     307        return kTRUE;
     308    }
     309
     310
     311   return kFALSE;
     312
     313}
     314
     315
     316
     317
     318
    239319
    240320
     
    340420    //    ls2: ls^2
    341421    //     ct: Cos(ZA.) - Cos(fThetaOffset)
    342     //    dd2: (DIST- fDistOffset)^2
     422    //    dd2: DIST^2 - fDistOffset^2
     423
     424  //     lconc: log10CONC   
     425
    343426
    344427    Double_t limit;
     
    488571
    489572
     573
     574    // The parametrization of upper cut in CONC is just provisional.
     575    // The aim is to get rid of triangular events.
     576    // The parametrization is
     577    // Log10(Conc)Uppercut = m*(log(Size) - log(SizeOffset)) + b
     578
     579
     580    Double_t log10conc = TMath::Log10(conc);
     581
     582
    490583    // in the parameterization of the cuts
    491584    // the dist parameter used is the one computed from the source position,
     
    508601    const Double_t asym    = asym0   * fMm2Deg;
    509602
    510 
     603   
     604    // MuuonLikeEvent is a variable which is set
     605    // to kFALSE/kTRUE by the filter cuts implemented
     606    // in MuonFilter function. If kTRUE, the event
     607    // will be rejected.
     608
     609// DM:
     610    //Bool_t MuonLikeEvent = MuonFilter (size, length);
     611   
    511612   
    512613    // computation of the cut limits
     
    524625    Double_t asymlow  = CtsMCut (fSuper->GetAsymLo(),   dmls, dmcza, dmls2, dd2);
    525626
    526     Double_t concup   = CtsMCut (fSuper->GetConcUp(),   dmls, dmcza, dmls2, dd2);
    527     Double_t conclow  = CtsMCut (fSuper->GetConcLo(),   dmls, dmcza, dmls2, dd2);
     627    Double_t log10concup   = CtsMCut (fSuper->GetConcUp(),   dmls, dmcza, dmls2, dd2);
     628    Double_t log10conclow  = CtsMCut (fSuper->GetConcLo(),   dmls, dmcza, dmls2, dd2);
    528629
    529630    Double_t  leakageup = CtsMCut (fSuper->GetLeakage1Up(),dmls, dmcza, dmls2, dd2);
     
    533634    Double_t lengthoverwidth = length/width;
    534635
     636
     637    // Cut in CONC is parametrized
     638
    535639    Double_t hadronness = 0.0;
    536640
    537641
    538     // If the cut values computed before for Dist, Length and Width exceed the upper limits set
    539     // by the user through function MSupercutsCalcONOFF::SetHillasDistLengthWidthUpperLowerLimits,
    540     // (or the default values defined in constructor); such cut values are replaced by the
     642    // If the cut values computed before for Dist,
     643    //  Length and Width exceed the upper limits set
     644    // by the user through function
     645    // MSupercutsCalcONOFF::SetHillasDistLengthWidthUpperLowerLimits,
     646    // (or the default values defined in constructor); such cut values
     647    // are replaced by the
    541648    // the limit values
    542649
     
    576683      }
    577684
     685    if (log10concup > fLog10ConcUpperLimit)
     686      {
     687        log10concup = fLog10ConcUpperLimit;
     688      }
     689
    578690   
    579691    if (leakageup > fLeakage1UpperLimit)
     
    582694      }
    583695   
    584 
    585    
     696   
    586697
    587698
    588699    // Upper cut in leakage 1 also implemented
    589700
    590    
    591    
    592    
    593701    if (//
    594702        newdist > distup ||
     
    600708        leakage > leakageup ||
    601709        leakage < leakagelow ||
    602         lengthoverwidth < fLengthOverWidthUpperLimit
     710        lengthoverwidth < fLengthOverWidthUpperLimit ||
    603711        //
    604712        //asym    < asymup &&
     
    608716        //dist    > distlow &&
    609717       
    610         //conc    < concup &&
    611         //conc    > conclow && 
     718        log10conc    > log10concup
     719        // log10conc    < log10conclow
    612720        //
     721        //log10conc    > log10concup ||
     722        // DM MuonLikeEvent // if KTRUE, event is rejected     
    613723        )
    614724         
     
    641751        ShowerParams[3] = newdist;
    642752        ShowerParams[4] = asym;
    643         ShowerParams[5] = conc;
     753        ShowerParams[5] = log10conc;
    644754        ShowerParams[6] = leakage;
    645755        ShowerParams[7] = theta;
     
    666776        SuperCutParams[6] = asymup;
    667777        SuperCutParams[7] = asymlow;
    668         SuperCutParams[8] = concup;
    669         SuperCutParams[9] = conclow;
     778        SuperCutParams[8] = log10concup;
     779        SuperCutParams[9] = log10conclow;
    670780        SuperCutParams[10] = leakageup;
    671781        SuperCutParams[11] = leakagelow;
  • trunk/MagicSoft/Mars/mtemp/mmpi/SupercutsONOFFClasses/MSupercutsCalcONOFF.h

    r4411 r6232  
    9494    Double_t fWidthUpperLimit;
    9595
     96    Double_t fLog10ConcUpperLimit;
     97
    9698     Double_t fDistLowerLimit;
    9799    Double_t fLengthLowerLimit;
     
    149151                                                  Double_t widthlow);
    150152
     153      // Function implementing a filter for muon induced images
     154// (Keichi cuts October 18th 2004)
     155
     156  Bool_t MuonFilter (Double_t Size,
     157                     Double_t Length);
     158
     159
     160
    151161    void InitMapping(MHMatrix *mat); // use quantity ThetaOrig.fVal at theta
    152162   
  • trunk/MagicSoft/Mars/mtemp/mmpi/SupercutsONOFFClasses/MTSupercutsApplied.cc

    r4411 r6232  
    130130    StringVector[6] = ("AsymUp/D");
    131131    StringVector[7] = ("AsymLow/D");
    132     StringVector[8] = ("ConcUp/D");
    133     StringVector[9] = ("ConcLow/D");
     132    StringVector[8] = ("Log10ConcUp/D");
     133    StringVector[9] = ("Log10ConcLow/D");
    134134    StringVector[10] = ("LeakageUp/D");
    135135    StringVector[11] = ("LeakageLow/D");
     
    155155    StringVector2[3] = ("NewDist/D");
    156156    StringVector2[4] = ("Asym/D");
    157     StringVector2[5] = ("Conc/D");
     157    StringVector2[5] = ("Log10Conc/D");
    158158    StringVector2[6] = ("Leakage/D");
    159159    StringVector2[7] = ("Theta/D");
  • trunk/MagicSoft/Mars/mtemp/mmpi/macros/SuperCutsONOFFMacro.C

    r4412 r6232  
    88gROOT -> Reset();
    99
    10 void SuperCutsONOFFMacro()
     10void SuperCutsONOFFMacroNew()
    1111{
    1212  gLog.SetNoColors();
     
    1515
    1616    // From magicserv01
    17     TString  ONDataFilename("/.magic/data16a/mazin/data/Mrk421/2004_04_22/4slices/Hillas_20040422_4sl_time_clean/Mrk421_*_HillasON.root");
     17//    TString  ONDataFilename("/.magic/data16a/mazin/data/mcdata/2004_06_30/*HillasON.root");
    1818                                       
    19     TString  OFFDataFilename("/.magic/data16a/mazin/data/Mrk421/2004_04_22/4slices/Hillas_20040422_4sl_time_clean/Mrk421_*_HillasOFF.root");
    20 
    21    
    22     TString PathForFiles ("/mnt/magicserv01/scratch/David/SillyTestForCommiting_July20_2004/");
    23 
    24    
    25     // **********************************************
     19//    TString  OFFDataFilename("/.magic/data16a/mazin/data/mcdata/2004_06_30/*HillasOFF.root");
     20
     21//   TString  ONDataFilename("/.magic/magicserv01/scratch/Period21/HillasFiles/CrabNebula_20040914.root");
     22//   TString  ONDataFilename("/.magic/magicserv01/scratch/Period21/HillasFiles/CrabNebula_20040921.root");
     23//   TString  ONDataFilename("/.magic/magicserv01/scratch/Period21/HillasFiles/CrabNebula_20040922.root");
     24
     25//   TString  ONDataFilename("/.magic/magicserv01/scratch/Period21/HillasFiles/CrabNebula_2004*.root");
     26 
     27   //TString  ONDataFilename("/.magic/data22a/mazin/HillasParam/Period24/2004_12_18/Crab_40_25_1_AbsHillas.root");
     28   TString  ONDataFilename("/.magic/data04a/nadia/Analysis/starfiles/DF/40_25_1/CrabNebula_*_Hillas_Abs.root");
     29   TString  OFFDataFilename("/.magic/data04a/nadia/Analysis/starfiles/DF/40_25_1/OffCrab1*Hillas_Abs.root");
     30   //TString ONDataFilename("/.magic/data03a/mazin/results/2004_09_21/CrabNebulaNadiaHillas.root"); 
     31   //TString OFFDataFilename("/.magic/data03a/mazin/results/2004_09_21/CrabNebulaNadiaHillas.root"); 
     32
     33
     34    TString PathForFiles ("/.magic/data03a/mazin/results/Crab/DF/CrabNadia/SuperCuts/Size2000/"); 
     35
     36//    TString PathForFiles ("/.magic/magicserv01/scratch/hbartko/SC/"); // was: Daniel/SuperCuts/Mrk421/2004_04_22/4slices_3520_nc/E800_1200_Opt_MC_Test/
     37
     38   
     39
    2640    // Boolean variables defining the job of the
    2741    // macro
    28     // **********************************************
    29 
    30    
    31 
     42 
    3243    // Boolean variable that decides wether data is read from files specified above
    3344    // (ON/OFF) or read from already existing Matrices (which are obviously stored
     
    3647    // values must be specified by user.
    3748
    38     // kTRUE reads alredy existing matrices, and kFALSE read data and produce matrices.
    39 
    40     Bool_t ReadMatrixFromRootFiles = kTRUE;
     49    Bool_t ReadMatrixFromRootFiles = kTRUE; 
    4150 
    4251
    43     // Boolean variable that controls wether to use the
    44     // TRAIN sample or not.
    45    
    46 
    47     Bool_t TrainParams   = kTRUE; 
    48 
    49 
    50 
    51     // Variable that allows the user to skip the optimization on the
    52     // train sample. If optimization is skipped (value kTRUE), the
    53     // previously optimized supercuts (stored in root file
    54     // which is called OptSCParametersONOFFThetaRangeXXXXXmRad.root, and located
    55     // in the directory specified by variable PathForFiles) are used
    56     // on the train and/or the test sample.
    57 
    58     // If value kFALSE, the cuts are optimized.
     52    // Boolean variable that controls the supercuts
     53    // optimization using the training sample
    5954    // The optimized cuts will be  written in root file
    6055    // located in directory specified before. Name of
    6156    // the root files is created automatically.
    62    
    63     Bool_t SkipOptimization = kTRUE;
    64  
    65 
    66 
    67 
     57
     58    Bool_t TrainParams   = kTRUE; 
     59
     60    // Variable that allows the user to skip the optimization on the
     61    // train sample. If optimization is skipped (value kTRUE), the
     62    // previously optimized supercuts (stored in root file) are used
     63    // on the train sample.
     64   
     65    Bool_t SkipOptimization = kFALSE;
    6866
    6967    // Boolean variable that allows the user to write the initial parameters
    7068    // into the root file that will be used to store the optimum cuts.
    71     // If ApplyInitialParams = kTRUE , the initial
    72     // parameters are written into this root file, and they
    73     // will be applied to the data (TRAIN/TEST )
    74     // IF NO OPTIMIZATION PROCEDURE IS PERFORMED.
    75    
    76     // If cuts are optimized (ie, variable SkipOptimization = kFALSE),
    77     // the cuts applied to the data are the optimized cuts.
     69    // If fUseInitialSCParams = kTRUE , parameters are written.
     70    // In this way, the initial SC parameters can be applied on the data (train/test)
     71   
     72    // The initial parameters are ONLY written to the root file if
     73    // there is NO SC params optimization, i.e., if variable
     74    // fSkipOptimization = kTRUE;
    7875 
    79     // NOTE: be aware that, if ApplyInitialSCParams = kTRUE and
    80     // there was a root file with the optimized cuts
     76    // NOTE: be aware that, if there was a file with the optimized cuts
    8177    // (previously computed), it will be overwritten with the initial
    82     // SC parameters.
    83 
    84     Bool_t ApplyInitialSCParams = kTRUE;
    85    
    86 
    87 
    88     // Boolean variable that controls wether to use the
    89     // TEST sample or not.
    90    
    91 
    92     Bool_t TestParams = kFALSE; 
     78    // SC parameters. This is something that I WILL HAVE TO CHANGE IN
     79    // future. Yet for the time being...
     80
     81    Bool_t UseInitialSCParams = kTRUE;  // was kFALSE
     82   
     83
     84
     85    // Variable that decides whether the optimized cuts are used
     86    // in the test sample.
     87
     88    Bool_t TestParams = kTRUE; 
    9389
    9490
     
    107103
    108104    // Fraction of ON events used for the training/testing
    109     Double_t whichfractiontrain = 0.999;
    110     Double_t whichfractiontest = 0.001;
     105    Double_t whichfractiontrain = 0.5;
     106    Double_t whichfractiontest = 0.5;
    111107
    112108    // Fraction of OFF events used for the training/testing
    113     Double_t whichfractiontrainOFF = 0.999;
    114     Double_t whichfractiontestOFF = 0.001;
     109    Double_t whichfractiontrainOFF = 0.5;
     110    Double_t whichfractiontestOFF = 0.5;
    115111
    116112
     
    126122    // Alpha value (degrees) below which signal is expected
    127123   
    128     Double_t alphasig = 12;
     124    Double_t alphasig = 9;
    129125
    130126    // Definition of alpha bkg region (where no signal is expected)
     
    135131    // Definition of the Size range
    136132
    137     Double_t SizeLow = 800;
    138     Double_t SizeUp = 1200;
    139 //    Double_t SizeUp = 1000000;
    140    
     133    Double_t SizeLow = 2000;
     134//    Double_t SizeUp =  2000;
     135    Double_t SizeUp = 1000000;
     136   
     137    Double_t LeakageMax = 0.04;
     138    Double_t DistMax    = 1.5;
     139    Double_t DistMin    = 0.1;
    141140
    142141     // Definition of binning of alpha plots
     
    167166    // in the computation of teh dynamical cuts that take place within
    168167    // class MCT1SupercutsCalc
    169     Bool_t NotUseTheta = kTRUE; // kTRUE renoves theta from the parameterization of cuts
     168    //Bool_t NotUseTheta = kFALSE; // kTRUE removes theta from the parameterization of cuts
     169    Bool_t NotUseTheta = kTRUE; // kTRUE removes theta from the parameterization of cuts
    170170
    171171    // Boolean variable used to decide wether to use dynamical cuts or static cuts
     
    192192    // Boolean variable used to decide wether initial parameters are
    193193    // read from ascii file or not. If kTRUE, parameters are retrieved
    194     // from ascii file. Otherwise, default parameters from MSupercuts
    195     // class are used.
     194    // from ascii file.
    196195
    197196    Bool_t ReadInitParamsFromAsciiFile = kTRUE;
     
    219218      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421DynCutsOnSize.txt"};
    220219      // {"../InitialSCParametersSteps/StartingValuesForOptimizationMkn421DynCutsOnSizeAndDist.txt"};
    221       {"mtemp/mmpi/asciifiles/OptimizedMkn421DynCutsGridWithSelected22pointsMay19.txt"};
     220      //{"mtemp/mmpi/asciifiles/OptimizedMkn421DynCutsGridWithSelected22pointsMay19.txt"};
     221      //{"/home/pcmagic16/mazin/mars/Mars260804/mtemp/mmpi/asciifiles/StartingValuesForSmallSizes1.txt"};
     222//      {"/home/pcmagic16/mazin/mars/Mars260804/mtemp/mmpi/asciifiles/SmallSizeRFStartValues.txt"};
     223      {"/home/pcmagic16/mazin/mars/Mars260804/mtemp/mmpi/asciifiles/OptimizedCrabFrom2000_a.txt"};
     224      //{"mtemp/mmpi/asciifiles/VeryGoodDummyValuesWithConcCut.txt"};
    222225
    223226
     
    244247
    245248    TArrayD CosThetaRangeVector(2);
    246     CosThetaRangeVector[0] = 0.0;
     249    CosThetaRangeVector[1] = 0.866;  //30
     250    //CosThetaRangeVector[1] = 0.74;  //42
     251    //CosThetaRangeVector[1] = 0.82;  //35
     252    //CosThetaRangeVector[1] = 0;
    247253    //CosThetaRangeVector[1] = 0.825;
    248254    //CosThetaRangeVector[2] = 0.921;
    249255    //CosThetaRangeVector[3] = 0.961;
    250     CosThetaRangeVector[1] = 1.0;
     256    CosThetaRangeVector[0] = 0.5;  // 60
     257    //CosThetaRangeVector[0] = 1;
    251258
    252259
     
    303310
    304311    FindSupercuts.SetSizeRange(SizeLow, SizeUp);
     312    FindSupercuts.SetFilters(LeakageMax, DistMax, DistMin);
     313//    FindSupercuts.SetFilters(0.001, 1.5, 0.2);
    305314
    306315
     
    318327    FindSupercuts.SetTrainParameters(TrainParams);
    319328    FindSupercuts.SetSkipOptimization(SkipOptimization);
    320     FindSupercuts.SetUseInitialSCParams(ApplyInitialSCParams);
     329    FindSupercuts.SetUseInitialSCParams(UseInitialSCParams);
    321330
    322331    FindSupercuts.SetTestParameters(TestParams);
Note: See TracChangeset for help on using the changeset viewer.