From d173b97ed2bd51520d1e2ea4abb42abf8363ce45 Mon Sep 17 00:00:00 2001
From: Alexander Voigt <Alexander.Voigt@desy.de>
Date: Thu, 21 Jul 2016 16:37:46 +0200
Subject: [PATCH] allow user to input CKM and PMNS matrices

---
 templates/librarylink.cpp.in | 38 +++++++++++++++++++++++++++++++++---
 templates/librarylink.m.in   | 20 +++++++++++++++++++
 2 files changed, 55 insertions(+), 3 deletions(-)

diff --git a/templates/librarylink.cpp.in b/templates/librarylink.cpp.in
index 4c2704d2f..a99f11962 100644
--- a/templates/librarylink.cpp.in
+++ b/templates/librarylink.cpp.in
@@ -295,11 +295,23 @@ void put_sm_input_parameters(const Model_data& data, MLINK link)
    MLPutRuleTo(link, data.qedqcd.displayMs2GeV(), "ms2GeV");
    MLPutRuleTo(link, data.qedqcd.displayMcMc(), "mcmc");
 
+   const flexiblesusy::CKM_parameters ckm(data.qedqcd.displayCKM());
+   MLPutRuleTo(link, ckm.theta_12, "CKMTheta12");
+   MLPutRuleTo(link, ckm.theta_13, "CKMTheta13");
+   MLPutRuleTo(link, ckm.theta_23, "CKMTheta23");
+   MLPutRuleTo(link, ckm.delta   , "CKMDelta");
+
+   const flexiblesusy::PMNS_parameters pmns(data.qedqcd.displayPMNS());
+   MLPutRuleTo(link, pmns.theta_12, "PMNSTheta12");
+   MLPutRuleTo(link, pmns.theta_13, "PMNSTheta13");
+   MLPutRuleTo(link, pmns.theta_23, "PMNSTheta23");
+   MLPutRuleTo(link, pmns.delta   , "PMNSDelta");
+   MLPutRuleTo(link, pmns.alpha_1 , "PMNSAlpha1");
+   MLPutRuleTo(link, pmns.alpha_2 , "PMNSAlpha2");
+
    MLPutRuleTo(link, data.physical_input.get(Physical_input::alpha_em_0), "alphaEm0");
    MLPutRuleTo(link, data.physical_input.get(Physical_input::mh_pole), "Mh");
 
-   // put CKM and PMNS matrices
-
    MLEndPacket(link);
 }
 
@@ -389,7 +401,7 @@ Model_data fill(double* pars, mint npars)
 {
    Model_data data;
 
-   const mint n_settings = 17, n_sm_parameters = 19, n_input_pars = @numberOfInputParameters@;
+   const mint n_settings = 17, n_sm_parameters = 29, n_input_pars = @numberOfInputParameters@;
    const mint n_total = n_settings + n_sm_parameters + n_input_pars;
 
    if (npars != n_total)
@@ -444,6 +456,26 @@ Model_data fill(double* pars, mint npars)
    data.qedqcd.setMass(softsusy::mCharm, pars[c]);
    data.qedqcd.setMcMc(pars[c++]);
 
+   {
+      flexiblesusy::CKM_parameters ckm;
+      ckm.theta_12 = pars[c++];
+      ckm.theta_13 = pars[c++];
+      ckm.theta_23 = pars[c++];
+      ckm.delta    = pars[c++];
+      data.qedqcd.setCKM(ckm);
+   }
+
+   {
+      flexiblesusy::PMNS_parameters pmns;
+      pmns.theta_12 = pars[c++];
+      pmns.theta_13 = pars[c++];
+      pmns.theta_23 = pars[c++];
+      pmns.delta    = pars[c++];
+      pmns.alpha_1  = pars[c++];
+      pmns.alpha_2  = pars[c++];
+      data.qedqcd.setPMNS(pmns);
+   }
+
    data.physical_input.set(Physical_input::alpha_em_0, pars[c++]);
    data.physical_input.set(Physical_input::mh_pole, pars[c++]);
 
diff --git a/templates/librarylink.m.in b/templates/librarylink.m.in
index 0407544be..c84feb35b 100644
--- a/templates/librarylink.m.in
+++ b/templates/librarylink.m.in
@@ -57,6 +57,16 @@ fsDefaultSMParameters = {
     mu2GeV -> 0.0024,       (* SMINPUTS[22] *)
     ms2GeV -> 0.104,        (* SMINPUTS[23] *)
     mcmc -> 1.27,           (* SMINPUTS[24] *)
+    CKMTheta12 -> 0,
+    CKMTheta13 -> 0,
+    CKMTheta23 -> 0,
+    CKMDelta -> 0,
+    PMNSTheta12 -> 0,
+    PMNSTheta13 -> 0,
+    PMNSTheta23 -> 0,
+    PMNSDelta -> 0,
+    PMNSAlpha1 -> 0,
+    PMNSAlpha2 -> 0,
     alphaEm0 -> 1/137.035999074,
     Mh -> 125.09
 };
@@ -114,6 +124,16 @@ FS@ModelName@OpenHandle[OptionsPattern[]] :=
             OptionValue[mu2GeV],
             OptionValue[ms2GeV],
             OptionValue[mcmc],
+            OptionValue[CKMTheta12],
+            OptionValue[CKMTheta13],
+            OptionValue[CKMTheta23],
+            OptionValue[CKMDelta],
+            OptionValue[PMNSTheta12],
+            OptionValue[PMNSTheta13],
+            OptionValue[PMNSTheta23],
+            OptionValue[PMNSDelta],
+            OptionValue[PMNSAlpha1],
+            OptionValue[PMNSAlpha2],
             OptionValue[alphaEm0],
             OptionValue[Mh]