Changeset 149 for GLBFrejus/HEAD/src


Ignore:
Timestamp:
May 30, 2006, 11:09:42 AM (18 years ago)
Author:
campagne
Message:

new version

Location:
GLBFrejus/HEAD/src
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • GLBFrejus/HEAD/src/discoCP.cxx

    r2 r149  
    3030  7/7/05  JEC introduce the Mass hierachy and Octant degeneracy from Th.Schwetz
    3131  8/7/05  JEC introduce data card management
     32  17/5/06 JEC add  TEST_OCTANT to avoid pb if theta23=pi/4.
    3233 */
    3334
    … …  
    3536#define DELTA_0   0.0
    3637#define DELTA_PI  M_PI
     38#define TEST_OCTANT false
    3739
    3840
    … …  
    4749//---------------------------------------------------------------------
    4850bool test_th23(glb_params fit_values) {
    49   double fit_th23 = glbGetOscParams(fit_values, GLB_THETA_23);
    50  
    51   // testing if fit_th23 is at the same side as the true value
    52  
    53   return  (
    54            (WRONG_TH23 == 0 &&
    55             (TH23 > M_PI/4. && fit_th23 > M_PI/4.) ||
    56             (TH23 < M_PI/4. && fit_th23 < M_PI/4.))
    57            ||
    58            (WRONG_TH23 == 1 &&
    59             (TH23 > M_PI/4. && fit_th23 < M_PI/4.) ||
    60             (TH23 < M_PI/4. && fit_th23 > M_PI/4.))
    61            );
     51  if(TEST_OCTANT) {
     52    double fit_th23 = glbGetOscParams(fit_values, GLB_THETA_23);
     53   
     54    // testing if fit_th23 is at the same side as the true value
     55   
     56    return  (
     57             (WRONG_TH23 == 0 &&
     58              (TH23 > M_PI/4. && fit_th23 > M_PI/4.) ||
     59              (TH23 < M_PI/4. && fit_th23 < M_PI/4.))
     60             ||
     61             (WRONG_TH23 == 1 &&
     62              (TH23 > M_PI/4. && fit_th23 < M_PI/4.) ||
     63              (TH23 < M_PI/4. && fit_th23 > M_PI/4.))
     64             );
     65  }
     66  return true;
    6267}//same_th23
    6368
    … …  
    101106
    102107  fileDataCard >> DUMMYSTR >> WRONG_TH23;
    103   std::cout << DUMMYSTR << WRONG_TH23 << std::endl;
     108  std::cout << DUMMYSTR << " " << WRONG_TH23 << std::endl;
    104109  fileDataCard >> DUMMYSTR >> WRONG_HIER;
    105   std::cout << DUMMYSTR << WRONG_HIER << std::endl;
     110  std::cout << DUMMYSTR << " " << WRONG_HIER << std::endl;
    106111  fileDataCard >> DUMMYSTR >> TH12; 
    107   std::cout << DUMMYSTR << TH12 << std::endl;
     112  std::cout << DUMMYSTR << " " << TH12 << std::endl;
    108113  fileDataCard >> DUMMYSTR >> TH23; 
    109   std::cout << DUMMYSTR << TH23 << std::endl;
     114  std::cout << DUMMYSTR << " " << TH23 << std::endl;
    110115  fileDataCard >> DUMMYSTR >> DMQ21; 
    111   std::cout << DUMMYSTR << DMQ21 << std::endl;
     116  std::cout << DUMMYSTR << " " << DMQ21 << std::endl;
    112117  fileDataCard >> DUMMYSTR >> DMQ31; 
    113   std::cout << DUMMYSTR << DMQ31 << std::endl;
     118  std::cout << DUMMYSTR << " " << DMQ31 << std::endl;
    114119  fileDataCard >> DUMMYSTR >> ROOTFILENAME; 
    115120
    … …  
    127132
    128133  fileDataCard >> DUMMYSTR >> ND; 
    129   std::cout << DUMMYSTR << ND  << std::endl;
     134  std::cout << DUMMYSTR << " "  << ND  << std::endl;
    130135  fileDataCard >> DUMMYSTR >> NTH; 
    131   std::cout << DUMMYSTR << NTH  << std::endl;
     136  std::cout << DUMMYSTR << " " << NTH  << std::endl;
    132137  fileDataCard >> DUMMYSTR >> LOG10MIN; 
    133   std::cout << DUMMYSTR << LOG10MIN << std::endl;
     138  std::cout << DUMMYSTR << " " << LOG10MIN << std::endl;
    134139  fileDataCard >> DUMMYSTR >> LOG10MAX; 
    135   std::cout << DUMMYSTR << LOG10MAX << std::endl;
     140  std::cout << DUMMYSTR << " "  << LOG10MAX << std::endl;
    136141
    137142
    … …  
    184189  // minimizer
    185190  glb_params input_errors = glbAllocParams();
     191  //JEC 17-05-06 use 10% except for Dm2_21 4%
     192//   glbDefineParams(input_errors,
     193//                TH12*0.05, 0., TH23*0.30, 0.,
     194//                DMQ21*0.05, fabs(DMQ31)*0.30);
    186195  glbDefineParams(input_errors,
    187                   TH12*0.05, 0., TH23*0.30, 0.,
    188                   DMQ21*0.05, fabs(DMQ31)*0.30);
     196                  TH12*0.10, 0., TH23*0.10, 0.,
     197                  DMQ21*0.04, fabs(DMQ31)*0.10);
    189198  glbSetDensityParams(input_errors, 0.05, GLB_ALL);
    190199  glbSetInputErrors(input_errors);
    … …  
    213222   
    214223    for(int j = 0; j < ND; j++) {
    215       //scan de 0 a +Pi pour delta_CP
    216       delta =  M_PI * double(j) / (ND-1);
     224      //scan de 0 a +2Pi pour delta_CP
     225      delta =  2.0 * M_PI * double(j) / (ND-1);
    217226
    218227      std::cout << "Consider (k,j) = (" << k << "," << j << ")" << std::endl;

  • GLBFrejus/HEAD/src/discoTheta13.cxx

    r4 r149  
    2424  experiment may be fit by Theta13 = 0
    2525  the true values are used to build the referenced experiments
    26   use chi2Theta to marginalized over all parameters except Theta
     26  use chi2Delta to marginalized over all parameters except delta
     27 
    2728 
    2829  11/7/05 JEC creation
    29   31/10/05 JEC correct the comment chi2Theta
    3030 */
    3131
    3232
    3333#define TH13_0   0.0
     34#define TEST_OCTANT false
    3435
    3536int WRONG_TH23, WRONG_HIER; //condition on octant/hierarchy
    … …  
    4344//---------------------------------------------------------------------
    4445bool test_th23(glb_params fit_values) {
    45   double fit_th23 = glbGetOscParams(fit_values, GLB_THETA_23);
    46  
    47   // testing if fit_th23 is at the same side as the true value
    48  
    49   return  (
    50            (WRONG_TH23 == 0 &&
    51             (TH23 > M_PI/4. && fit_th23 > M_PI/4.) ||
    52             (TH23 < M_PI/4. && fit_th23 < M_PI/4.))
    53            ||
    54            (WRONG_TH23 == 1 &&
    55             (TH23 > M_PI/4. && fit_th23 < M_PI/4.) ||
    56             (TH23 < M_PI/4. && fit_th23 > M_PI/4.))
    57            );
     46  if(TEST_OCTANT) {
     47    double fit_th23 = glbGetOscParams(fit_values, GLB_THETA_23);
     48   
     49    // testing if fit_th23 is at the same side as the true value
     50   
     51    return  (
     52             (WRONG_TH23 == 0 &&
     53              (TH23 > M_PI/4. && fit_th23 > M_PI/4.) ||
     54              (TH23 < M_PI/4. && fit_th23 < M_PI/4.))
     55             ||
     56             (WRONG_TH23 == 1 &&
     57              (TH23 > M_PI/4. && fit_th23 < M_PI/4.) ||
     58              (TH23 < M_PI/4. && fit_th23 > M_PI/4.))
     59             );
     60  }
     61  return true;
    5862}//same_th23
    5963

  • GLBFrejus/HEAD/src/sensi.cxx

    r2 r149  
    2525  use glbChiDeltaTheta to marginalized other all parameters except of
    2626  the current delta and theta
     27
     28  24/10/05 (JEC)
     29  . transfert the #defined variables to Data Cards
     30  . polish the code
     31
    2732 */
    28 //FIXME: transfert the #defined variables to Data Cards
    29 #define TH12 ( 0.5 * asin(sqrt(0.82)) )
    30 #define TH23 ( M_PI/4.0 )
    31 
    32 #define TH13   0.0 
    33 #define DELTA  0.0
    34 
    35 #define DMQ21 8.2e-5
    36 #define DMQ31 2.5e-3
    37 
    38 
    39 int WRONG_TH23, WRONG_HIER;
    40 
     33
     34float TH12, TH23,DMQ21,DMQ31; //theta12,theta23,Dm^2_21,Dm^2_31
     35
     36//The true value Theta13 and Delta_CP
     37#define TH13_0   0.0 
     38#define DELTA_0  0.0
     39
     40int ND; //number of bins to scan delta parameter [-pi,+pi]
     41int NTH; //number of bins to scan theta13 parameter
     42float LOG10MIN; //minimal Log10(sin^2(2theta13))
     43float LOG10MAX; //minimal Log10(sin^2(2theta13))
    4144
    4245/**************************************/
    … …  
    4649int main(int argc, char *argv[])
    4750{
     51  //Process data card file
     52  std::string ROOTDIR = getenv("GLBFREJUSROOT");
     53  std::string DATACARD;
     54  DATACARD = ROOTDIR + "/run/sensi.data";
     55  std::ifstream fileDataCard;
     56  std::cout << "Open DataCard file : " << DATACARD <<std::endl;
     57  fileDataCard.open(DATACARD.data());
     58  std::string DUMMYSTR; //comment string of the variables
     59
     60  std::string ROOTFILENAME;
     61
     62  fileDataCard >> DUMMYSTR >> TH12; 
     63  std::cout << DUMMYSTR << TH12 << std::endl;
     64  fileDataCard >> DUMMYSTR >> TH23; 
     65  std::cout << DUMMYSTR << TH23 << std::endl;
     66  fileDataCard >> DUMMYSTR >> DMQ21; 
     67  std::cout << DUMMYSTR << DMQ21 << std::endl;
     68  fileDataCard >> DUMMYSTR >> DMQ31; 
     69  std::cout << DUMMYSTR << DMQ31 << std::endl;
     70  fileDataCard >> DUMMYSTR >> ROOTFILENAME; 
     71  ROOTFILENAME += ".root";
     72  std::cout << DUMMYSTR <<ROOTFILENAME  << std::endl;
     73
     74  fileDataCard >> DUMMYSTR >> ND; 
     75  std::cout << DUMMYSTR << ND  << std::endl;
     76  fileDataCard >> DUMMYSTR >> NTH; 
     77  std::cout << DUMMYSTR << NTH  << std::endl;
     78  fileDataCard >> DUMMYSTR >> LOG10MIN; 
     79  std::cout << DUMMYSTR << LOG10MIN << std::endl;
     80  fileDataCard >> DUMMYSTR >> LOG10MAX; 
     81  std::cout << DUMMYSTR << LOG10MAX << std::endl;
     82
     83
    4884  //init AIDA for Histo/Tuple
    4985  AIDA::IAnalysisFactory* aida = AIDA_createAnalysisFactory();
    … …  
    5793  std::string opts = "export=root";
    5894  AIDA::ITree* fTree =
    59     treeFactory->create("SPLDeltaTheta2+8OldFluxCorr.root","root",false,true,opts);
     95    treeFactory->create(ROOTFILENAME,"root",false,true,opts);
    6096  delete treeFactory;
    6197
    … …  
    88124  // true values (reference for future Chi2 computations)
    89125  glb_params true_values = glbAllocParams();
    90   glbDefineParams(true_values, TH12, TH13, TH23, DELTA, DMQ21, DMQ31);
     126  glbDefineParams(true_values, TH12, TH13_0, TH23, DELTA_0, DMQ21, DMQ31);
    91127  glbPrintParams(stdout,true_values);
    92128 
    … …  
    96132  // set Starting Values
    97133  glb_params start_values = glbAllocParams();
    98   glbDefineParams(start_values, TH12, TH13, TH23, DELTA, DMQ21, DMQ31);
     134  glbDefineParams(start_values, TH12, TH13_0, TH23, DELTA_0, DMQ21, DMQ31);
    99135  glbSetStartingValues(start_values);
    100136  glbPrintParams(stdout,start_values);
    … …  
    103139  // init test values
    104140  glb_params test_values = glbAllocParams();
    105   glbDefineParams(test_values, TH12, TH13, TH23, DELTA, DMQ21, DMQ31);
     141  glbDefineParams(test_values, TH12, TH13_0, TH23, DELTA_0, DMQ21, DMQ31);
    106142  glbPrintParams(stdout,test_values);
    107143
    … …  
    122158
    123159
    124 #define ND  31 
    125 #define NTH 41
    126 #define LOG10MIN -4.0
    127 #define LOG10MAX -1.0
    128 
    129160  std::cout << "Start to loop...." << std::endl;
    130161
    … …  
    147178      glbSetOscParams(test_values, th, GLB_THETA_13);
    148179     
     180      std::cout << "Consider (k,j) = (" << k << "," << j << ")" << std::endl;
    149181      //Param correlation
    150182      res = glbChiThetaDelta(test_values,fit_values,0);
    151 
    152       std::cout << "delta j:" << j << " theta13 k:" << k
    153                 << " chi2 : " << res
    154                 << std::endl;
    155 
     183     
     184
     185      //For sensitivity plot store the outcome of the fit
    156186      fit_theta12 = glbGetOscParams(fit_values, GLB_THETA_12);
    157187      fit_theta13 = glbGetOscParams(fit_values, GLB_THETA_13);
Note: See TracChangeset for help on using the changeset viewer.