loon: NC::OscProb::FourFlavorBase Class Reference

00228 {
00229   if(intType == NCType::kNC){
00230     switch(mode){
00231     case NCType::kNuMuToNuMu:
00232       return 1-TransitionProbability(NCType::kNuMuToNuS, intType, baseline, trueE);
00233     case NCType::kNuMuToNuTau:
00234     case NCType::kNuMuToNuE:
00235       return 0;
00236     case NCType::kNuEToNuE:
00237       return 1;
00238     case NCType::kNuMuToNuS:
00239       // Fall through to actually calculate the probability
00240       break;
00241     default:
00242       assert(0 && "Unknown mode");
00243     }
00244   }
00245 
00246   fLoverE = baseline/trueE;
00247 
00248   // Precalculate sines and cosines
00249   const double s13 = TMath::Sin(Theta13());
00250   const double s23 = TMath::Sin(Theta23());
00251   const double s14 = TMath::Sin(Theta14());
00252   const double s24 = TMath::Sin(Theta24());
00253   const double s34 = TMath::Sin(Theta34());
00254 
00255   const double c13 = TMath::Cos(Theta13());
00256   const double c23 = TMath::Cos(Theta23());
00257   const double c14 = TMath::Cos(Theta14());
00258   const double c24 = TMath::Cos(Theta24());
00259   const double c34 = TMath::Cos(Theta34());
00260 
00261   // Precalculate complex phases
00262   const TComplex emd1(1, -Delta1(), true); // 'true' selects polar coordinates
00263   const TComplex epd1 = TComplex::Conjugate(emd1);
00264   const TComplex emd2(1, -Delta2(), true);
00265   const TComplex epd2 = TComplex::Conjugate(emd2);
00266 
00267   // Matrix elements
00268   const TComplex ue3 = emd1*c14*s13;
00269 
00270   const TComplex ue4 = s14;
00271 
00272   const TComplex umu3 = emd1*emd2*-s13*s14*s24 + c13*s23*c24;
00273 
00274   const TComplex umu4 = emd2*c14*s24;
00275 
00276   const TComplex utau3_1 = emd1*-s13*s14*c24*s34;
00277   const TComplex utau3_2 = epd2*-c13*s23*s24*s34;
00278   const TComplex utau3 = c13*c23*c34 + utau3_1 + utau3_2;
00279 
00280   const TComplex utau4 = c14*c24*s34;
00281 
00282   const TComplex us3_1 = emd1*-s13*s14*c24*c34;
00283   const TComplex us3_2 = epd2*-c13*s23*s24*c34;
00284   const TComplex us3 = -c13*c23*s34 + us3_1 + us3_2;
00285 
00286   const TComplex us4 = c14*c24*c34;
00287 
00288   // Find the square of the moduli for the elements
00289   const double ue3Sqr   = ue3.Rho2();
00290   const double ue4Sqr   = ue4.Rho2();
00291   const double umu3Sqr  = umu3.Rho2();
00292   const double umu4Sqr  = umu4.Rho2();
00293   const double utau3Sqr = utau3.Rho2();
00294   const double utau4Sqr = utau4.Rho2();
00295   const double us3Sqr   = us3.Rho2();
00296   const double us4Sqr   = us4.Rho2();
00297 
00298   // Check unitarity
00299   const double epsilon = 1e-5;
00300   const double col3Sum  = ue3Sqr+umu3Sqr+utau3Sqr+us3Sqr;
00301   const double col4Sum = ue4Sqr+umu4Sqr+utau4Sqr+us4Sqr;
00302   if(TMath::Abs(col3Sum-1) > epsilon ||
00303      TMath::Abs(col4Sum-1) > epsilon)
00304     MSG("NCOscProb", Msg::kWarning) << "Warning: unitarity violated!!!!"
00305                                     << endl;
00306 
00307   const TComplex umu4conj = TComplex::Conjugate(umu4);
00308 
00309   // Calculate some cross-terms
00310   const TComplex umuue   = umu4conj*ue4  *umu3*TComplex::Conjugate(ue3);
00311   const TComplex umuutau = umu4conj*utau4*umu3*TComplex::Conjugate(utau3);
00312   const TComplex umuus   = umu4conj*us4  *umu3*TComplex::Conjugate(us3);
00313 
00314   double p  = 0;
00315 
00316 #ifdef CASEMUTOALPHA
00317 #error Someone already defined CASEMUTOALPHA, have to choose a different name
00318 #endif
00319 
00320 
00321 #define CASEMUTOALPHA(alpha) \
00322     p  = umu3Sqr*u##alpha##3Sqr*sinSqrDelta31(); \
00323     p += umu4Sqr*u##alpha##4Sqr*sinSqrDelta41(); \
00324     p += (umuu##alpha).Re()*(sinSqrDelta31() - sinSqrDelta43() + sinSqrDelta41()); \
00325     p *= 4; \
00326     p += 2.*(umuu##alpha).Im()*(sin2Delta31() - sin2Delta41() + sin2Delta43()); \
00327     return p
00328 
00329   switch(mode){
00330   case NCType::kNuMuToNuMu:
00331     p  = umu3Sqr*(1-umu3Sqr-umu4Sqr)*sinSqrDelta31();
00332     p += umu4Sqr*(1-umu3Sqr-umu4Sqr)*sinSqrDelta41();
00333     p += umu4Sqr*umu3Sqr*sinSqrDelta43();
00334     p *= 4.;
00335     return 1 - p;
00336 
00337   case NCType::kNuMuToNuTau:
00338     CASEMUTOALPHA(tau);
00339 
00340   case NCType::kNuMuToNuE:
00341     CASEMUTOALPHA(e);
00342 
00343   case NCType::kNuMuToNuS:
00344     CASEMUTOALPHA(s);
00345 
00346   case NCType::kNuEToNuE:
00347     // OK since there are very few beam nues anyway?
00348     return 1;
00349 
00350   default:
00351     assert(0 && "Unknown mode");
00352   }
00353 
00354 #undef CASEMUTOALPHA
00355 }


Public Member Functions
	FourFlavorBase (NCType::EOscModel mod)
virtual	~FourFlavorBase ()
double	TransitionProbability (NCType::EOscMode mode, NCType::EEventType interactionType, double baseline, double trueE) const
Protected Member Functions
virtual double	Theta13 () const =0
	$\theta_{13}$
virtual double	Theta23 () const =0
	$\theta_{23}$
virtual double	Theta14 () const =0
	$\theta_{14}$
virtual double	Theta24 () const =0
	$\theta_{24}$
virtual double	Theta34 () const =0
	$\theta_{34}$
virtual double	Delta1 () const =0
	$\delta_1$
virtual double	Delta2 () const =0
	$\delta_2$
virtual double	sinSqrDelta31 () const =0
	$\sin^2\Delta_{31}$
virtual double	sin2Delta31 () const =0
	$\sin2\Delta_{31}$
virtual double	sinSqrDelta41 () const =0
	$\sin^2\Delta_{41}$
virtual double	sin2Delta41 () const =0
	$\sin2\Delta_{41}$
virtual double	sinSqrDelta43 () const =0
	$\sin^2\Delta_{43}$
virtual double	sin2Delta43 () const =0
	$\sin2\Delta_{43}$
Protected Attributes
double	fLoverE

NC::OscProb::FourFlavorBase Class Reference

Public Member Functions

Protected Member Functions

Protected Attributes

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation