Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Class List | File List | Namespace Members | Class Members | File Members | Related Pages

neugen_inputs.cxx

Go to the documentation of this file.
00001 
00013 #include <fstream>
00014 #include <sstream>
00015 #include <iomanip>
00016 
00017 #include "neugen_inputs.h"
00018 
00019 ClassImp(neugen_inputs)
00020 
00021 using std::endl;
00022 using std::setw;
00023 using std::ios;
00024 using std::setiosflags;
00025 using std::setfill;
00026 using std::ofstream;
00027 using std::ostringstream;
00028 
00029 //____________________________________________________________________________
00030 ostream & operator << (ostream & stream, const neugen_inputs & conf)
00031 {
00032   conf.print(stream);
00033   return stream;
00034 }
00035 //____________________________________________________________________________
00036 neugen_inputs::neugen_inputs()
00037 {
00038   init();
00039 }
00040 //____________________________________________________________________________
00041 neugen_inputs::neugen_inputs(const neugen_inputs * inputs)
00042 {
00043   init();
00044 
00045    // neugen card params
00046 
00047   _nbins             = inputs->_nbins;       
00048   _xsec_type_code    = inputs->_xsec_type_code;
00049   _emin              = inputs->_emin;
00050   _emax              = inputs->_emax;
00051   _e                 = inputs->_e;
00052   _plot_var_code     = inputs->_plot_var_code;
00053   _flux_id_code      = inputs->_flux_id_code;
00054   _plot_range_code   = inputs->_plot_range_code;
00055   _plot_var_min      = inputs->_plot_var_min;
00056   _plot_var_max      = inputs->_plot_var_max;
00057   _nu_type_code      = inputs->_nu_type_code;
00058   _wcurrent_code     = inputs->_wcurrent_code;
00059   _target_code       = inputs->_target_code;
00060   _final_state_code  = inputs->_final_state_code;
00061   _init_state_code   = inputs->_init_state_code;
00062   _process_mask_code = inputs->_process_mask_code;
00063   _cut_var_code      = inputs->_cut_var_code;
00064   _cut_var_min       = inputs->_cut_var_min;
00065   _cut_var_max       = inputs->_cut_var_max;
00066   _qel_sum           = inputs->_qel_sum;
00067   _res_sum           = inputs->_res_sum;
00068   _dis_sum           = inputs->_dis_sum;
00069 
00070   _qel_bit_in_mask   = inputs->_qel_bit_in_mask;
00071   _res_bit_in_mask   = inputs->_res_bit_in_mask;
00072   _dis_bit_in_mask   = inputs->_dis_bit_in_mask;
00073 
00074    // aux
00075 
00076   _fin_p             = inputs->_fin_p;
00077   _fin_n             = inputs->_fin_n;       
00078   _fin_pi_plus       = inputs->_fin_pi_plus;   
00079   _fin_pi_0          = inputs->_fin_pi_0;   
00080   _fin_pi_minus      = inputs->_fin_pi_minus;  
00081   _xsec_type_str     = inputs->_xsec_type_str;
00082   _plot_var_str      = inputs->_plot_var_str;
00083   _flux_id_str       = inputs->_flux_id_str;
00084   _plot_range_str    = inputs->_plot_range_str;
00085   _nu_type_str       = inputs->_nu_type_str;
00086   _wcurrent_str      = inputs->_wcurrent_str;
00087   _target_str        = inputs->_target_str;
00088   _final_state_str   = inputs->_final_state_str;
00089   _init_state_str    = inputs->_init_state_str;
00090   _cut_var_str       = inputs->_cut_var_str;
00091 }
00092 //____________________________________________________________________________
00093 neugen_inputs::~neugen_inputs()
00094 {
00095 
00096 }
00097 //____________________________________________________________________________
00098 void neugen_inputs::write_neugen_input_card(const char * filename) const
00099 {
00100 // This methods writes out the object state in the form of 'data cards' that
00101 // NeuGEN understands.
00102 //  
00103   ofstream data_card(filename);
00104 
00105   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00106             << _nbins
00107             << " \\\\ nbins:  number of points in plot "
00108             << endl;
00109   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00110             << _xsec_type_code
00111             << " \\\\ xsec type: 1=total, 2=differential "
00112             << endl;
00113   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00114             << _emin
00115             << " \\\\ xsec type=1 - Emin "
00116             << endl;
00117   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00118             << _emax
00119             << " \\\\             - Emax "
00120             << endl;
00121   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00122             << _e
00123             << " \\\\ xsec type=2 - E "
00124             << endl;
00125   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00126             << _plot_var_code
00127             << " \\\\ plot variable "
00128             << endl;
00129   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00130             << _flux_id_code
00131             << " \\\\ flux id (1=ANL, 2=GGM, 3=BNL, 4=BEBC) "
00132             << endl;;
00133   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00134             << _plot_range_code
00135             << " \\\\ plot  range (1=auto, 2=custom) "
00136             << endl;
00137   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00138             << _plot_var_min
00139             << " \\\\ plot range - min "
00140             << endl;
00141   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00142             << _plot_var_max
00143             << " \\\\            - max "
00144             << endl;
00145   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00146             << _nu_type_code
00147             << " \\\\ neutrino (nue/bar=5/6, numu/bar=7/8, nutau/bar=9/10 "
00148             << endl;
00149   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00150             << _wcurrent_code
00151             << " \\\\ weak current (1=CC, 2=NC, 3=BOTH) "
00152             << endl;
00153   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00154             << _target_code
00155             << " \\\\ target: nucleus / particle code, -1 for isoscalar "
00156             << endl;
00157   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00158             << _final_state_code
00159             << " \\\\ final state - in the form pn+-0 "
00160             << endl;
00161   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00162             << _init_state_code
00163             << " \\\\ initial state - (1=v-p, 2=v-n, 3=vbar-p, 4=vbar-n "
00164             << endl;
00165   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00166             << _process_mask_code
00167             << " \\\\ process mask: bits for qel, res, dis, coh "
00168             << endl;
00169   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00170             << _cut_var_code
00171             << " \\\\ cuts variable (0=none, 1=|q^2|, 2=W, 3=x, 4=y "
00172             << endl;
00173   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00174             << _cut_var_min
00175             << " \\\\ cut variable - min "
00176             << endl;
00177   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00178             << _cut_var_max
00179             << " \\\\             - max "
00180             << endl;
00181   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00182             << _qel_sum 
00183             << " \\\\ qelsum: >0 means add qel channel "
00184             << endl;
00185   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00186             << _res_sum 
00187             << " \\\\ ressum: >0 means add all res channels "
00188             << endl;
00189   data_card << setiosflags(ios::left) << setfill(' ') << setw(9)
00190             << _dis_sum 
00191             << " \\\\ dissum: >0 means add all dis "
00192             << endl;
00193 }
00194 //____________________________________________________________________________
00195 interaction neugen_inputs::get_interaction(void) const
00196 {
00197   flavor_t      f = flavor::get_from_code(_nu_type_code);
00198   nucleus_t     n = e_free;
00199   ccnc_t        c = ccnc::get_from_code(_wcurrent_code);
00200   init_state_t  i = init_state::get_init_state_from_code(_init_state_code);
00201 
00202   interaction inter(f, n, c, i);
00203 
00204   return inter;
00205 }
00206 //____________________________________________________________________________
00207 final_state neugen_inputs::get_final_state(void) const
00208 {
00209   final_state state;
00210 
00211   state.setFinalState(_fin_p, _fin_n, _fin_pi_plus, _fin_pi_minus, _fin_pi_0);
00212 
00213   return state;
00214 }
00215 //____________________________________________________________________________
00216 neugen_cuts neugen_inputs::get_cuts(void) const
00217 {
00218   kinematic_variable_t kvid = 
00219                   kinematic_variable::get_kin_var_id_from_code(_cut_var_code);
00220 
00221   bool sumQel = (_qel_sum == 1);
00222   bool sumRes = (_res_sum == 1);
00223   bool sumDis = (_dis_sum == 1);
00224 
00225   neugen_cuts cuts(kvid, _cut_var_min, _cut_var_max, 
00226                                   _process_mask_code, sumQel, sumRes, sumDis);
00227 
00228   return cuts;
00229 }
00230 //____________________________________________________________________________
00231 void neugen_inputs::set_nbins(int nbins)
00232 {
00233   _nbins = nbins;
00234 }
00235 //____________________________________________________________________________
00236 void neugen_inputs::set_qel_bit_in_mask(bool on)
00237 {
00238   ( on ) ? _qel_bit_in_mask = 1 : _qel_bit_in_mask = 0;
00239 
00240   compute_process_mask();
00241 }
00242 //____________________________________________________________________________
00243 void neugen_inputs::set_res_bit_in_mask(bool on)
00244 {
00245   ( on ) ? _res_bit_in_mask = 1 : _res_bit_in_mask = 0;
00246 
00247   compute_process_mask();
00248 }
00249 //____________________________________________________________________________
00250 void neugen_inputs::set_dis_bit_in_mask(bool on)
00251 {
00252   ( on ) ? _dis_bit_in_mask = 1 : _dis_bit_in_mask = 0;
00253 
00254   compute_process_mask();
00255 }
00256 //____________________________________________________________________________
00257 void neugen_inputs::set_xsec_type(string xsec_type)
00258 {
00259   _xsec_type_str  = xsec_type;
00260   _xsec_type_code = neugen_xsec_type_code(xsec_type);
00261 }
00262 //____________________________________________________________________________
00263 void neugen_inputs::set_e_min(float e_min)
00264 {
00265   _emin = e_min;
00266 }
00267 //____________________________________________________________________________
00268 void neugen_inputs::set_e_max(float e_max)
00269 {
00270   _emax = e_max;
00271 }
00272 //____________________________________________________________________________
00273 void neugen_inputs::set_e( float e)
00274 {
00275   _e = e;
00276 }
00277 //____________________________________________________________________________
00278 void neugen_inputs::set_plot_var(string plot_variable)
00279 {
00280   _plot_var_str  = plot_variable;
00281   _plot_var_code = neugen_variable_code(plot_variable);
00282 }
00283 //____________________________________________________________________________
00284 void neugen_inputs::set_flux(string flux)
00285 {
00286   _flux_id_str  = flux;
00287   _flux_id_code = neugen_flux_code(flux);
00288 }
00289 //____________________________________________________________________________
00290 void neugen_inputs::set_range(string range)
00291 {
00292   _plot_range_str  = range;
00293   _plot_range_code = neugen_plot_range_code(range);
00294 }
00295 //____________________________________________________________________________
00296 void neugen_inputs::set_plot_var_min(float var_min)
00297 {
00298   _plot_var_min = var_min;
00299 }
00300 //____________________________________________________________________________
00301 void neugen_inputs::set_plot_var_max(float var_max)
00302 {
00303   _plot_var_max = var_max;
00304 }
00305 //____________________________________________________________________________
00306 void neugen_inputs::set_neutrino(string neutrino)
00307 {
00308   _nu_type_str  = neutrino;
00309   _nu_type_code = neugen_neutrino_code(neutrino);
00310 }
00311 //____________________________________________________________________________
00312 void neugen_inputs::set_wcurrent(string wcurrent)
00313 {
00314   _wcurrent_str  = wcurrent;
00315   _wcurrent_code = neugen_weak_current_code(wcurrent);
00316 }
00317 //____________________________________________________________________________
00318 void neugen_inputs::set_target(string /*target*/)
00319 {
00320   _target_str  = ""; // unused
00321   _target_code = 0;  // unused
00322 }
00323 //____________________________________________________________________________
00324 void neugen_inputs::set_cut_var(string cut_variable)
00325 {
00326   _cut_var_str  = cut_variable;
00327   _cut_var_code = neugen_variable_code(cut_variable);
00328 }
00329 //____________________________________________________________________________
00330 void neugen_inputs::set_cut_var_min(float var_min)
00331 {
00332   _cut_var_min = var_min;
00333 }
00334 //____________________________________________________________________________
00335 void neugen_inputs::set_cut_var_max(float var_max)
00336 {
00337   _cut_var_max = var_max;
00338 }
00339 //____________________________________________________________________________
00340 void neugen_inputs::set_inclusive(bool on)
00341 {
00342   _inclusive = on;
00343   
00344   set_qel_sum(on);
00345   set_res_sum(on);
00346   set_dis_sum(on);
00347 }
00348 //____________________________________________________________________________
00349 void neugen_inputs::set_qel_sum(bool on)
00350 {
00351   _qel_sum  = bool2int(on);
00352 }
00353 //____________________________________________________________________________
00354 void neugen_inputs::set_res_sum(bool on)
00355 {
00356   _res_sum  = bool2int(on);
00357 }
00358 //____________________________________________________________________________
00359 void neugen_inputs::set_dis_sum(bool on)
00360 {
00361   _dis_sum  = bool2int(on);
00362 }
00363 //____________________________________________________________________________
00364 void neugen_inputs::set_final_state(string fin_state)
00365 {
00366   _final_state_str  = fin_state;
00367   _final_state_code = neugen_final_state_code(fin_state);
00368 }
00369 //____________________________________________________________________________
00370 void neugen_inputs::set_initial_state(string init_state)
00371 {
00372   _init_state_str  = init_state;
00373   _init_state_code = neugen_initial_state_code(init_state);
00374 }  
00375 //____________________________________________________________________________
00376 int neugen_inputs::neugen_xsec_type_code(string xsec_type)
00377 {
00378   if (xsec_type.find("differential") != string::npos) return 2;
00379   else return 1;
00380 }
00381 //____________________________________________________________________________
00382 int neugen_inputs::neugen_flux_code(string flux)
00383 {
00384   if      (flux.find("ANL")  != string::npos) return 1;
00385   else if (flux.find("GGM")  != string::npos) return 2;
00386   else if (flux.find("BNL")  != string::npos) return 3;
00387   else if (flux.find("BEBC") != string::npos) return 4;
00388   else return 0;  
00389 }
00390 //____________________________________________________________________________
00391 int neugen_inputs::neugen_plot_range_code(string range)
00392 {
00393   if (range.find("custom") != string::npos) return 2;
00394   else return 1;
00395 }
00396 //____________________________________________________________________________
00397 int neugen_inputs::neugen_neutrino_code(string neutrino)
00398 {
00399   if      (neutrino.find("nu_e")       != string::npos) return  5;
00400   else if (neutrino.find("nu_e_bar")   != string::npos) return  6;
00401   else if (neutrino.find("nu_mu")      != string::npos) return  7;
00402   else if (neutrino.find("nu_mu_bar")  != string::npos) return  8;
00403   else if (neutrino.find("nu_tau")     != string::npos) return  9;
00404   else if (neutrino.find("nu_tau_bar") != string::npos) return 10;
00405   else                                                  return  0;
00406 }
00407 //____________________________________________________________________________
00408 int neugen_inputs::neugen_weak_current_code(string wcurrent)
00409 {
00410   if (wcurrent.find("+") != string::npos)       return 3;
00411   else if (wcurrent.find("NC") != string::npos) return 2;
00412   else                                          return 1;
00413 }
00414 //____________________________________________________________________________
00415 string neugen_inputs::neugen_final_state_code(string fin_state)
00416 {
00417 // build neugen's final state - in the form (pn+-0)
00418 
00419   _fin_p        = 0;   // init
00420   _fin_n        = 0;
00421   _fin_pi_plus  = 0;
00422   _fin_pi_0     = 0;
00423   _fin_pi_minus = 0;
00424 
00425   _fin_p        = n_matches(fin_state, "p ");
00426   _fin_n        = n_matches(fin_state, "n ");
00427   _fin_pi_plus  = n_matches(fin_state, "pi(+)");
00428   _fin_pi_0     = n_matches(fin_state, "pi(0)");
00429   _fin_pi_minus = n_matches(fin_state, "pi(-)");
00430 
00431   ostringstream code;
00432 
00433   code << _fin_p << _fin_n << _fin_pi_plus << _fin_pi_minus << _fin_pi_0;
00434 
00435   return code.str();
00436 }
00437 //____________________________________________________________________________
00438 int neugen_inputs::neugen_initial_state_code(string init_state)
00439 {
00440   if      (init_state.find("nu + p")     != string::npos) return 1;
00441   else if (init_state.find("nu + n")     != string::npos) return 2;
00442   else if (init_state.find("nu_bar + p") != string::npos) return 3;
00443   else if (init_state.find("nu_bar + n") != string::npos) return 4;
00444   else if (init_state.find("nu + N")     != string::npos) return 5; // ????
00445   else if (init_state.find("nu_bar + N") != string::npos) return 6; // ????
00446   else                                              return 0;
00447 }
00448 //____________________________________________________________________________
00449 int neugen_inputs::neugen_variable_code(string cut_variable)
00450 {
00451   if      (cut_variable.find("none")  != string::npos) return 0;
00452   else if (cut_variable.find("|q^2|") != string::npos) return 1;
00453   else if (cut_variable.find("W")     != string::npos) return 2;
00454   else if (cut_variable.find("x")     != string::npos) return 3;
00455   else if (cut_variable.find("y")     != string::npos) return 4;
00456   else                                                 return 0;
00457 }      
00458 //____________________________________________________________________________
00459 void neugen_inputs::compute_process_mask(void)
00460 {
00461   int qel, dis, res;
00462   
00463   (_qel_bit_in_mask == 1) ? qel = 0 : qel = 1;
00464   (_res_bit_in_mask == 1) ? res = 0 : res = 1;
00465   (_dis_bit_in_mask == 1) ? dis = 0 : dis = 1;
00466   
00467   _process_mask_code = qel + 2 * res + 4 * dis;
00468 }
00469 //____________________________________________________________________________
00470 int neugen_inputs::n_matches(string input, string pattern)
00471 {
00472   // if max = 1 then this is ok - do something more generic later
00473   
00474   if (input.find(pattern) != string::npos) return 1;
00475   else return 0;
00476 
00477   /*
00478   string::size_type pos = 0;
00479   int n=0;
00480   
00481   while( (pos = input.find_first_of(pattern, pos)) != string::npos ) {
00482 
00483       n++;
00484       input.erase(pos, pattern.length());
00485   }
00486   return n;
00487   */
00488 }
00489 //____________________________________________________________________________
00490 int neugen_inputs::bool2int(bool on)
00491 {
00492   if (on) return 1;
00493   else    return 0;
00494 }
00495 //____________________________________________________________________________
00496 void neugen_inputs::init(void)
00497 {
00498   //-- init neugen cards variables
00499   
00500   _nbins             = 0;
00501   _xsec_type_code    = 0;
00502   _emin              = 0;
00503   _emax              = 0;
00504   _e                 = 0;
00505   _plot_var_code     = 0;
00506   _flux_id_code      = 0;
00507   _plot_range_code   = 0;
00508   _plot_var_min      = 0;
00509   _plot_var_max      = 0;
00510   _nu_type_code      = 0;
00511   _wcurrent_code     = 0;
00512   _target_code       = 0;
00513   _final_state_code  = "00000";
00514   _init_state_code   = 0;
00515   _process_mask_code = 0;
00516   _cut_var_code      = 0;
00517   _cut_var_min       = 0;
00518   _cut_var_max       = 0;
00519   _qel_sum           = 0;
00520   _res_sum           = 0;
00521   _dis_sum           = 0;
00522 
00523   _qel_bit_in_mask   = 0;
00524   _res_bit_in_mask   = 0;
00525   _dis_bit_in_mask   = 0;
00526 
00527   //-- init auxiliary variables
00528 
00529   _fin_p             = 0;
00530   _fin_n             = 0;
00531   _fin_pi_plus       = 0;
00532   _fin_pi_0          = 0;
00533   _fin_pi_minus      = 0;
00534   
00535   _xsec_type_str     = "";
00536   _plot_var_str      = "";
00537   _flux_id_str       = "";
00538   _plot_range_str    = "";
00539   _nu_type_str       = "";
00540   _wcurrent_str      = "";
00541   _target_str        = "";
00542   _final_state_str   = "";
00543   _init_state_str    = "";
00544   _cut_var_str       = "";    
00545 }
00546 //____________________________________________________________________________
00547 void neugen_inputs::print(ostream & stream) const
00548 {
00549   stream << "number of bins =  " << _nbins             << endl;
00550   stream << "xsec type =       " << _xsec_type_code    << endl;
00551   stream << "E min =           " << _emin              << endl;
00552   stream << "E max =           " << _emax              << endl;
00553   stream << "E =               " << _e                 << endl;
00554   stream << "plot var =        " << _plot_var_code     << endl;
00555   stream << "flux id =         " << _flux_id_code      << endl;
00556   stream << "plot range =      " << _plot_range_code   << endl;
00557   stream << "plot var - min =  " << _plot_var_min      << endl;
00558   stream << "plot var - max =  " << _plot_var_max      << endl;
00559   stream << "neutrino type =   " << _nu_type_code      << endl;
00560   stream << "weak current =    " << _wcurrent_code     << endl;
00561   stream << "target =          " << _target_code       << endl;
00562   stream << "final state =     " << _final_state_code  << endl;
00563   stream << "initial state =   " << _init_state_code   << endl;
00564   stream << "process mask =    " << _process_mask_code << endl;
00565   stream << "qel bit in mask = " << _qel_bit_in_mask   << endl;
00566   stream << "res bit in mask = " << _res_bit_in_mask   << endl;
00567   stream << "dis bit in mask = " << _dis_bit_in_mask   << endl;
00568   stream << "cut variable =    " << _cut_var_code      << endl;
00569   stream << "cut var - min =   " << _cut_var_min       << endl;
00570   stream << "cut var - max =   " << _cut_var_max       << endl;
00571   stream << "qel sum =         " << _qel_sum           << endl;
00572   stream << "res sum =         " << _res_sum           << endl;
00573   stream << "dis sum =         " << _dis_sum           << endl;
00574 }
00575 //____________________________________________________________________________

Generated on Mon Mar 16 22:58:57 2009 for loon by doxygen 1.3.5