[342] | 1 | |
---|
[304] | 2 | #include "softwareGenerator.h" |
---|
| 3 | #include "abstractElement.h" |
---|
[312] | 4 | #include "generatorParticle.h" |
---|
| 5 | #include "mathematicalConstants.h" |
---|
| 6 | #include "PhysicalConstants.h" |
---|
[359] | 7 | #include "dataManager.h" |
---|
[304] | 8 | |
---|
| 9 | softwareGenerator::softwareGenerator() : abstractSoftware() |
---|
[386] | 10 | { |
---|
| 11 | nameOfSoftware_ = nomDeLogiciel("generator"); |
---|
| 12 | } |
---|
[304] | 13 | |
---|
[316] | 14 | softwareGenerator::softwareGenerator(string inputFileName, globalParameters* globals, dataManager* dt) : abstractSoftware(inputFileName, globals, dt) |
---|
[386] | 15 | { |
---|
| 16 | nameOfSoftware_ = nomDeLogiciel("generator"); |
---|
[418] | 17 | registerElement(nomdElements::beam,TBoolOk); |
---|
[386] | 18 | } |
---|
[304] | 19 | |
---|
[316] | 20 | bool softwareGenerator::createInputFile(particleBeam* beamBefore, unsigned int numeroDeb, unsigned int numeroFin, string workingDir) |
---|
[304] | 21 | { |
---|
[386] | 22 | if ( !initComputationLimits(numeroDeb,numeroFin) ) return false; |
---|
[304] | 23 | |
---|
[386] | 24 | if ( numeroDeb_ != numeroFin_ ) return false; |
---|
| 25 | |
---|
[304] | 26 | abstractElement* elPtr; |
---|
[386] | 27 | elPtr = dataManager_->getElementPointerFromNumero(numeroDeb_); |
---|
[304] | 28 | |
---|
[418] | 29 | if ( elPtr->getNomdElement().getElementType() != nomdElements::RFgun ) { |
---|
[313] | 30 | dataManager_->consoleMessage(" softwareGenerator::createInputFile : the element must be rfgun " ); |
---|
[304] | 31 | cerr << " softwareGenerator::createInputFile : the element must be rfgun" << endl; |
---|
| 32 | return false; |
---|
| 33 | } |
---|
[342] | 34 | |
---|
[304] | 35 | ofstream outfile; |
---|
[316] | 36 | string name = workingDir + inputFileName_; |
---|
[304] | 37 | outfile.open(name.c_str(), ios::out); |
---|
| 38 | if (!outfile) { |
---|
[313] | 39 | dataManager_->consoleMessage(" softwareGenerator::createInputFile : error opening output stream " ); |
---|
[304] | 40 | cerr << " error opening output stream " << name << endl; |
---|
| 41 | return false; |
---|
| 42 | } |
---|
| 43 | outfile << "&INPUT" << endl; |
---|
| 44 | string fichier = "'" + workingDir + "faisceau.ini" + "'"; |
---|
| 45 | outfile << " FNAME = " << fichier << endl; |
---|
| 46 | outfile << " Add=FALSE, N_add=0" << endl; |
---|
| 47 | outfile << " Species='electrons'" << endl; |
---|
[418] | 48 | // outfile << elPtr->generatorOutputFlow() << endl; |
---|
| 49 | |
---|
| 50 | outfile << elementsData(elPtr->parametersToSoftware()) << endl; |
---|
| 51 | |
---|
[304] | 52 | outfile << "/" << endl; |
---|
| 53 | outfile.close(); |
---|
[313] | 54 | dataManager_->consoleMessage("fichier input termine pour GENERATOR"); |
---|
[304] | 55 | return true; |
---|
| 56 | } |
---|
[312] | 57 | |
---|
[386] | 58 | bool softwareGenerator::execute(string workingDir) { |
---|
[342] | 59 | |
---|
| 60 | bool ExecuteStatus = true; |
---|
| 61 | |
---|
[312] | 62 | ostringstream sortie; |
---|
[386] | 63 | sortie << " EXECUTION DE GENERATOR DE l'ELEMENT " << numeroDeb_ << " A L'ELEMENT " << numeroFin_ << endl; |
---|
[312] | 64 | |
---|
| 65 | string generatorJob = workingDir + "generator"; |
---|
| 66 | generatorJob += string(" "); |
---|
[316] | 67 | generatorJob += workingDir + inputFileName_; |
---|
[312] | 68 | |
---|
| 69 | string resultOfRun; |
---|
[342] | 70 | bool success = launchJob(generatorJob,resultOfRun); |
---|
[312] | 71 | sortie << resultOfRun << endl; |
---|
[342] | 72 | if ( !success ) { |
---|
[312] | 73 | sortie << " launching of generator failed " << endl; |
---|
| 74 | ExecuteStatus = false; |
---|
[342] | 75 | } else { |
---|
[312] | 76 | cout << " execution generator MARCHE " << endl; |
---|
| 77 | sortie << resultOfRun; |
---|
[342] | 78 | string nameOut = workingDir + "generator.output"; |
---|
| 79 | ofstream outfile; |
---|
[312] | 80 | outfile.open(nameOut.c_str(), ios::out); |
---|
[342] | 81 | if (!outfile) { |
---|
| 82 | sortie << " error first opening transport output stream " << nameOut << endl; |
---|
| 83 | ExecuteStatus = false; |
---|
| 84 | } else { |
---|
[312] | 85 | // on copie la sortie dans un fichier 'generator.out' |
---|
| 86 | outfile << resultOfRun << endl; |
---|
| 87 | outfile.close(); |
---|
| 88 | } |
---|
| 89 | } |
---|
[342] | 90 | |
---|
| 91 | dataManager_->consoleMessage(sortie.str()); |
---|
[312] | 92 | return ExecuteStatus; |
---|
| 93 | } |
---|
| 94 | |
---|
[386] | 95 | bool softwareGenerator::buildBeamAfterElements( string workingDir) { |
---|
[312] | 96 | |
---|
| 97 | bool result = true; |
---|
[386] | 98 | if ( !ComputationLimitsOk() ) return false; |
---|
| 99 | |
---|
| 100 | if ( numeroDeb_ != numeroFin_ ) { |
---|
[313] | 101 | dataManager_->consoleMessage(" softwareGenerator::buildBeamAfterElements : only one element (rfgun) must be calculated " ); |
---|
| 102 | return false; |
---|
| 103 | } |
---|
[386] | 104 | // on initialise une nouvelle sortie diagnostic |
---|
| 105 | particleBeam* newDiag = dataManager_->updateCurrentDiagnostic(true); |
---|
[312] | 106 | |
---|
[386] | 107 | // beams.push_back(particleBeam()); |
---|
[312] | 108 | vector<double> centroid; |
---|
| 109 | bareParticle refPart; |
---|
| 110 | vector<bareParticle> particles; |
---|
| 111 | vector<bareParticle> particlesPassives; // on ne fait rien de ces particules pour l'instant |
---|
| 112 | if (beamFromGenerator(string("faisceau.ini"),workingDir, centroid, refPart,particles, particlesPassives )) { |
---|
[386] | 113 | newDiag->setWithParticles(centroid, refPart,particles); |
---|
[312] | 114 | } else { |
---|
[313] | 115 | dataManager_->consoleMessage(" softwareGenerator::buildBeamAfterElements : error " ); |
---|
[312] | 116 | result = false; |
---|
| 117 | } |
---|
| 118 | return result; |
---|
| 119 | } |
---|
| 120 | |
---|
| 121 | bool softwareGenerator::beamFromGenerator(string beamFileName, string workingDir, vector<double>& centroid, bareParticle& refPart,vector<bareParticle>& particles, vector<bareParticle>& passiveParticles ) { |
---|
| 122 | unsigned k; |
---|
| 123 | FILE* filefais; |
---|
| 124 | string nomfilefais = workingDir + beamFileName; |
---|
| 125 | cout << " nom fichier sortie generator : " << nomfilefais << endl; |
---|
| 126 | filefais = fopen(nomfilefais.c_str(), "r"); |
---|
| 127 | if ( filefais == (FILE*)0 ) { |
---|
[313] | 128 | dataManager_->consoleMessage(" softwareGenerator::beamFromGenerator : error openig the file " + nomfilefais ); |
---|
[312] | 129 | cerr << " beamFromGenerator() erreur a l'ouverture du fichier" << nomfilefais << endl;; |
---|
| 130 | return false; |
---|
| 131 | } |
---|
| 132 | else cout << " beamFromGenerator() : ouverture du fichier " << nomfilefais << endl; |
---|
| 133 | |
---|
| 134 | generatorParticle partic; |
---|
| 135 | std::vector<generatorParticle> faisceau; |
---|
| 136 | int nbProbPart =0; |
---|
| 137 | double timeRef = 0.0; |
---|
| 138 | double betagammaZRef = 0.0; |
---|
| 139 | // lecture part. de reference |
---|
| 140 | if ( partic.readFromGeneratorFile(filefais) > 0 ) |
---|
| 141 | { |
---|
| 142 | if ( partic.index != 1 ) { |
---|
[313] | 143 | dataManager_->consoleMessage(" softwareGenerator::beamFromGenerator : particles are not electrons, we have to reconsider this method beamFromGenerator" ); |
---|
| 144 | cout << " ERROR softwareGenerator::beamFromGenerator : particles are not electrons, we have to reconsider this method " << endl; |
---|
[312] | 145 | return false; |
---|
| 146 | } |
---|
| 147 | |
---|
| 148 | if ( partic.flag != -1 ) { |
---|
| 149 | cout << " ATTENTION softwareGenerator::beamFromGenerator : flag different de -1 " << endl; |
---|
| 150 | } |
---|
| 151 | if (fabs(partic.xx) > EPSILON || fabs(partic.yy) > EPSILON || fabs(partic.px) > EPSILON || fabs(partic.py) > EPSILON) { |
---|
| 152 | printf(" ATTENTION softwareGenerator::beamFromGenerator : part. reference douteuse \n"); |
---|
| 153 | partic.imprim(); |
---|
| 154 | } |
---|
| 155 | timeRef = partic.clock; |
---|
| 156 | TRIDVECTOR posRef(100.*partic.xx,100.*partic.yy,100.*partic.zz); // en cm |
---|
[356] | 157 | betagammaZRef = partic.pz/EREST_eV; |
---|
[312] | 158 | // l'impulsion donnee par generator est en eV/c |
---|
[356] | 159 | TRIDVECTOR betagammaRef(partic.px/EREST_eV , partic.py/EREST_eV, betagammaZRef); |
---|
[312] | 160 | refPart = bareParticle(posRef, betagammaRef); |
---|
| 161 | |
---|
| 162 | // seule la part. de reference a un pz absolu (les autres pat. on un pz relatif a la ref) |
---|
| 163 | // je mets ici a zero, pour homogeneiser |
---|
| 164 | partic.pz = 0.0; |
---|
| 165 | faisceau.push_back(partic); |
---|
| 166 | } |
---|
| 167 | |
---|
| 168 | while( partic.readFromGeneratorFile(filefais) > 0 ) { |
---|
| 169 | faisceau.push_back(partic); |
---|
| 170 | if ( partic.flag != -1 ) nbProbPart++; |
---|
| 171 | } |
---|
| 172 | |
---|
| 173 | if ( faisceau.size() == 0) |
---|
| 174 | { |
---|
[313] | 175 | dataManager_->consoleMessage(" softwareGenerator::beamFromGenerator : error no particle found" ); |
---|
[312] | 176 | cerr << " softwareGenerator::beamFromGenerator echec lecture " << endl; |
---|
| 177 | return false; |
---|
| 178 | } |
---|
| 179 | // pour l'instant on choisit un centroid nul; |
---|
| 180 | centroid.clear(); |
---|
| 181 | centroid = vector<double>(6,0.0); |
---|
| 182 | |
---|
| 183 | particles.clear(); |
---|
| 184 | passiveParticles.clear(); |
---|
| 185 | // particles.resize(faisceau.size() - nbProbPart, bareParticle()); |
---|
| 186 | // passiveParticles.resize(nbProbPart, bareParticle()); |
---|
| 187 | double x,y; |
---|
[417] | 188 | // double deltaz; |
---|
| 189 | double cdeltat; |
---|
[312] | 190 | TRIDVECTOR pos; |
---|
| 191 | TRIDVECTOR betagamma; |
---|
| 192 | double pxPart; |
---|
| 193 | double pyPart; |
---|
| 194 | double pzPartRel; |
---|
| 195 | double deltaPzPart; |
---|
| 196 | |
---|
| 197 | for ( k=0; k < faisceau.size(); k++) { |
---|
| 198 | |
---|
| 199 | pxPart = faisceau.at(k).px; |
---|
| 200 | pyPart = faisceau.at(k).py; |
---|
| 201 | pzPartRel = faisceau.at(k).pz; |
---|
| 202 | deltaPzPart = faisceau.at(k).pz; |
---|
| 203 | x=faisceau.at(k).xx; |
---|
| 204 | y=faisceau.at(k).yy; |
---|
| 205 | |
---|
| 206 | // tout ce qui suit sera a clarifier |
---|
[356] | 207 | double betaGammax = pxPart/EREST_eV; |
---|
| 208 | double betaGammay = pyPart/EREST_eV; |
---|
| 209 | double betaGammaz = betagammaZRef + pzPartRel/EREST_eV; |
---|
[312] | 210 | betagamma.setComponents(betaGammax, betaGammay, betaGammaz); |
---|
| 211 | double gamma = sqrt(1.0 + betagamma.norm2()); |
---|
| 212 | |
---|
[417] | 213 | // decalage temporel par rapport a la reference ? |
---|
| 214 | // on prend le faisceau sous la forme "a un instant donne" |
---|
[312] | 215 | |
---|
[417] | 216 | // deltat = faisceau.at(k).clock - timeRef; // nanoseondes |
---|
| 217 | cdeltat = CLIGHT_m_per_ns * (faisceau.at(k).clock - timeRef); // metres |
---|
| 218 | // double ds = CLIGHT_m_per_ns * deltat /gamma; |
---|
| 219 | // x += betaGammax * ds; // en metres |
---|
| 220 | // y += betaGammay * ds; |
---|
| 221 | |
---|
[312] | 222 | // ici on neglige la difference entre gamma de la part. et gamma de la ref. |
---|
[417] | 223 | // deltaz = (pzPartRel/EREST_eV) * CLIGHT_m_per_ns * deltat; // en metres |
---|
[312] | 224 | |
---|
[417] | 225 | // pos.setComponents(100.*x,100.*y,100.*deltaz); // en cm |
---|
| 226 | pos.setComponents(100.*x,100.*y,100.*cdeltat); // en cm |
---|
[312] | 227 | if ( faisceau.at(k).flag == -1 ) { |
---|
| 228 | particles.push_back(bareParticle(pos,betagamma)); |
---|
| 229 | } else { |
---|
| 230 | passiveParticles.push_back(bareParticle(pos,betagamma)); |
---|
| 231 | } |
---|
| 232 | } |
---|
| 233 | |
---|
| 234 | |
---|
| 235 | return true; |
---|
| 236 | } |
---|
[418] | 237 | |
---|
| 238 | |
---|
| 239 | string softwareGenerator::elementsData(const vector< pair<string, vector<string> > >& donnees) const { |
---|
| 240 | unsigned k; |
---|
| 241 | cout << " PASSAGE softwareGenerator::elementsData " << endl; |
---|
| 242 | if ( donnees.at(0).first != "labelsGenericSpecific" ) { |
---|
| 243 | cout << " softwareGenerator::elementsData ERROR : element badly defined " << endl; |
---|
| 244 | return string(); |
---|
| 245 | } |
---|
| 246 | string genericName = donnees.at(0).second.at(0); |
---|
| 247 | if ( genericName == "rfgun" ) return rfgunData(donnees); |
---|
| 248 | return string(); |
---|
| 249 | } |
---|
| 250 | |
---|
| 251 | string softwareGenerator::rfgunData(const vector< pair<string, vector<string> > >& donnees) const { |
---|
| 252 | |
---|
| 253 | cout << " PASSAGE softwareGenerator::rfgunData " << endl; |
---|
| 254 | string nmacrop = "0"; |
---|
| 255 | string sigma_t = "0.0"; |
---|
| 256 | string sigma_r = "0.0"; |
---|
| 257 | string emit_x = "0.0"; |
---|
| 258 | string emit_y = "0.0"; |
---|
| 259 | string E_cin = "0.0"; |
---|
| 260 | string sigma_E = "0.0"; |
---|
| 261 | string phaseStep = "0.0"; |
---|
| 262 | string totalCharge = "0.0"; |
---|
| 263 | |
---|
| 264 | unsigned k; |
---|
| 265 | for ( k=1; k < donnees.size(); k++) { |
---|
| 266 | if ( donnees.at(k).first == "nbMacroparticles" ) { |
---|
| 267 | nmacrop = donnees.at(k).second.at(0); |
---|
| 268 | } else if ( donnees.at(k).first == "sigmasTR" ) { |
---|
| 269 | sigma_t = donnees.at(k).second.at(0); |
---|
| 270 | sigma_r = donnees.at(k).second.at(1); |
---|
| 271 | } else if ( donnees.at(k).first == "emittancesXY" ) { |
---|
| 272 | emit_x = donnees.at(k).second.at(0); |
---|
| 273 | emit_y = donnees.at(k).second.at(1); |
---|
| 274 | } else if ( donnees.at(k).first == "kineticE" ) { |
---|
| 275 | E_cin = donnees.at(k).second.at(0); |
---|
| 276 | sigma_E = donnees.at(k).second.at(1); |
---|
| 277 | } else if ( donnees.at(k).first == "phaseStep" ) { |
---|
| 278 | phaseStep = donnees.at(k).second.at(0); |
---|
| 279 | } else if ( donnees.at(k).first == "totalCharge" ) { |
---|
| 280 | totalCharge = donnees.at(k).second.at(0); |
---|
| 281 | } |
---|
| 282 | } |
---|
| 283 | ostringstream sortie; |
---|
| 284 | |
---|
| 285 | |
---|
| 286 | sortie << "Ipart=" << nmacrop << endl; |
---|
| 287 | |
---|
| 288 | sortie << "Probe=.True." << endl; |
---|
| 289 | sortie << "Noise_reduc=.T" << endl; |
---|
| 290 | sortie << "Cathode=.T." << endl; |
---|
| 291 | sortie << "Q_total=" << totalCharge << endl; |
---|
| 292 | sortie << "Ref_zpos=0.0" << endl; |
---|
| 293 | sortie << "Ref_clock=0.0" << endl; |
---|
| 294 | sortie << "Ref_Ekin=" << E_cin << endl; // tjs en MeV |
---|
| 295 | sortie << "Dist_z='g'" << endl; |
---|
| 296 | sortie << "sig_clock=" << 1.0e-3*atof(sigma_t.c_str()) << endl; // passage de ps en ns |
---|
| 297 | sortie << "Dist_pz='g', sig_Ekin=" << 1000.*atof(sigma_E.c_str()) << ", emit_z=0.0 , cor_Ekin=0.0 " << endl; // passage en keV |
---|
| 298 | sortie << "Dist_x='gauss', sig_x=" << 10.*atof(sigma_r.c_str()) << endl; // passage en mm |
---|
| 299 | sortie << "Dist_px='g', Nemit_x=" << emit_x << ", cor_px=0.0" << endl; |
---|
| 300 | sortie << "Dist_y='gauss', sig_y=" << 10.*atof(sigma_r.c_str()) << endl; // passage en mm |
---|
| 301 | sortie << "Dist_py='g', Nemit_y=" << emit_y << ", cor_py=0.0" << endl; |
---|
| 302 | |
---|
| 303 | return sortie.str(); |
---|
| 304 | } |
---|