#include "ReactorFortran.hh"
#include "ReactorConstants.hh"
#include "ReactorEvent.hh"
#include "ReactorDetector.hh"
#include "ReactorXC.hh"
#include "ReactorNtuple.hh"
#include <iostream.h>
#include <fstream.h>
#include <iomanip.h>
#include <math.h>
#include <string.h>
#include <cstdio>
Go to the source code of this file.
Functions | |
int | main () |
Top level function of ReactorFsim. More... |
ReactorFsim simulates any number of spherical, three-zone liquid scintillator neutrino detectors.
Main author: Tim Bolton
Contributing authors: Erin Abouzaid, Matthew Worcester
Definition in file ReactorFsim.cpp.
|
Top level function of ReactorFsim. First sets the number of events, initializes paw, and defines each detector. The event loop creates an event and runs it through each detector added to the loop, so that each detector "sees" the exact same event each time. Finally, paw is closed and the hbook file created. Definition at line 29 of file ReactorFsim.cpp. References ReactorNtuple::Fill, ReactorDetector::LightsOut, ReactorDetector::SetParam_FastNeutronOption, ReactorDetector::SetParam_R0, ReactorDetector::SetParam_R1, and ReactorDetector::SetParam_R2.
00029 { 00030 00031 // ReactorNeutrons Neutrons; 00032 // 00033 // Set number of events, paw file name, and ntuple name 00034 // 00035 // defaults 00036 int nevents = 100; 00037 double r0=260; 00038 double r1=260; 00039 double rmaxgen=270; 00040 double r2=350; 00041 int fast=0; 00042 char* fileName="test.ntuple"; 00043 char* ntName="event data"; 00044 // 00045 // setup the event 00046 fstream fin("src/event_setup.txt", ios::in); 00047 //cout << "setup file: src/event_setup.txt" << endl; 00048 00049 char paramLine[256]; 00050 std::string dummy; 00051 bool done = false; 00052 std::string fastneutron; 00053 00054 while(!done){ 00055 00056 // '#' specifies a comment line -- ignore it 00057 if(fin.peek() == '#'){ 00058 00059 fin.getline(paramLine, 256); 00060 continue; 00061 } 00062 00063 // 'E'ND or 'e'nd specifies end of file -- set done flag to true 00064 if(fin.peek() == 'e' || fin.peek() == 'E'){ 00065 00066 done = true; 00067 continue; 00068 } 00069 00070 // 'N'EVENT or 'n'event specifies the number of events 00071 if(fin.peek() == 'n' || fin.peek() == 'N'){ 00072 00073 fin >> dummy; 00074 while(fin.peek() == ' ' || fin.peek() == '=') 00075 fin.ignore(1); 00076 00077 fin >> nevents; 00078 00079 if(nevents <= 0){ 00080 cout << " SETUP ERROR: number of events <= 0!" << endl; 00081 nevents = 1; 00082 } 00083 00084 continue; 00085 } 00086 00087 // 'O'UTPUTFILE or 'o'utputfile specifies the output filename 00088 if(fin.peek() == 'o' || fin.peek() == 'O'){ 00089 00090 fin >> dummy; 00091 while(fin.peek() == ' ' || fin.peek() == '=') 00092 fin.ignore(1); 00093 00094 char newFile[256]; 00095 fin >> newFile; 00096 fileName = newFile; 00097 00098 continue; 00099 } 00100 00101 // 'R'ADII or 'r'adii specifies the three detector radii: R0 R1 R2 00102 if(fin.peek() == 'r' || fin.peek() == 'R'){ 00103 00104 fin >> dummy; 00105 while(fin.peek() == ' ' || fin.peek() == '=') 00106 fin.ignore(1); 00107 00108 fin >> r0 >> r1 >> r2 >> rmaxgen; 00109 00110 continue; 00111 } 00112 00113 // 'F'ASTNEUTRON or 'f'astneutron turns on/off fast neutrons 00114 if(fin.peek() == 'f' || fin.peek() == 'F'){ 00115 00116 fin >> dummy; 00117 while(fin.peek() == ' ' || fin.peek() == '=') 00118 fin.ignore(1); 00119 00120 fin >> fastneutron; 00121 00122 if(fastneutron == "off"){ 00123 fast = 0; 00124 } 00125 else{ 00126 fast = 1; 00127 } 00128 00129 continue; 00130 } 00131 00132 // skip T, P. and D 00133 if(fin.peek() == 'T' || fin.peek() == 't' || 00134 fin.peek() == 'P' || fin.peek() == 'p' || 00135 fin.peek() == 'D' || fin.peek() == 'd'){ 00136 00137 fin.getline(paramLine, 256); 00138 continue; 00139 } 00140 00141 // ignore extra whitespace 00142 if(fin.peek() == ' ' || fin.peek() == '\n'){ 00143 00144 fin.ignore(1); 00145 continue; 00146 } 00147 } 00148 // 00149 cout << "ReactorFsim v0_0_2 by T. Bolton" << endl; 00150 cout << "===============================" << endl; 00151 cout << " " << nevents << " events" << endl; 00152 cout << " Output file: " << fileName << endl; 00153 cout << " R0 = " << r0 << ", R1 = " << r1 << ", R2 = " << r2 00154 << " (Rmaxgen = " << rmaxgen << ") cm" << endl; 00155 cout << " Fast neutrons are " << fastneutron << endl; 00156 cout << "===============================" << endl; 00157 // 00158 // init random number seed 00159 // 00160 int luxury = 2; 00161 int seed = 8796896; 00162 int k1=0; 00163 int k2=0; 00164 rluxgo_(luxury,seed,k1,k2); 00165 // 00166 // Paw setup. 00167 // 00168 hbook_init__(); 00169 // 00170 int iostat; 00171 int ntLun=99; 00172 int ntBuf=1024; 00173 char* dirName = "TEST"; 00174 char* chopt = "N"; 00175 hropen_(ntLun,dirName,fileName,chopt,ntBuf,iostat, 00176 strlen(dirName),strlen(fileName),strlen(chopt)); 00177 // 00178 // Set up simple detector models 00179 // Set up simple Ntuple and histogranms 00180 // See ReactorNtuple.hh for ntuple/histogram details. 00181 // 00182 ReactorDetector DetectorA; 00183 DetectorA.SetParam_FastNeutronOption(fast); 00184 DetectorA.SetParam_R0(r0); 00185 DetectorA.SetParam_R1(r1); 00186 DetectorA.SetParam_R2(r2); 00187 ReactorNtuple ModelA(1,"model A",dirName); 00188 // 00189 // event loop 00190 // 00191 for (int events = 0;events<nevents;events++){ 00192 // 00193 // create neutrino event 00194 // 00195 ReactorEvent Event(rmaxgen); 00196 // 00197 // create fast detector simulation 00198 // 00199 DetectorA.LightsOut(Event); 00200 ModelA.Fill(Event,DetectorA); 00201 } 00202 // 00203 // close up paw 00204 // 00205 hldir_(" "," ",strlen(" "),strlen(" ")); 00206 char dirName2[strlen(dirName)+2]; 00207 *dirName2='\0'; 00208 strcat(dirName2,"//"); 00209 strcat(dirName2,dirName); 00210 hcdir_(dirName2," ",strlen(dirName2),strlen(" ")); 00211 int outopt = 0; 00212 int icycle; 00213 hrout_(outopt,icycle," ",strlen(" ")); 00214 hrend_(dirName,strlen(dirName)); 00215 } |