File indexing completed on 2024-08-13 05:00:13
0001 #include "FWCore/Framework/interface/EventSetup.h"
0002 #include "FWCore/Framework/interface/ESHandle.h"
0003 #include "Geometry/Records/interface/PTrackerAdditionalParametersPerDetRcd.h"
0004 #include "Alignment/CommonAlignment/interface/SurveyDet.h"
0005 #include "Alignment/TrackerAlignment/interface/AlignableTracker.h"
0006 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeomBuilderFromGeometricDet.h"
0007
0008 #include "FWCore/Framework/interface/EventSetup.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011
0012 #include "DataFormats/DetId/interface/DetId.h"
0013 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0014 #include "Alignment/SurveyAnalysis/plugins/CreateSurveyRcds.h"
0015 #include "Geometry/GeometryAligner/interface/GeometryAligner.h"
0016 #include "CLHEP/Random/RandGauss.h"
0017
0018 CreateSurveyRcds::CreateSurveyRcds(const edm::ParameterSet& cfg)
0019 : tTopoToken_(esConsumes()),
0020 geomDetToken_(esConsumes()),
0021 ptpToken_(esConsumes()),
0022 ptitpToken_(esConsumes()),
0023 aliToken_(esConsumes()),
0024 aliErrToken_(esConsumes()) {
0025 m_inputGeom = cfg.getUntrackedParameter<std::string>("inputGeom");
0026 m_inputSimpleMis = cfg.getUntrackedParameter<double>("simpleMis");
0027 m_generatedRandom = cfg.getUntrackedParameter<bool>("generatedRandom");
0028 m_generatedSimple = cfg.getUntrackedParameter<bool>("generatedSimple");
0029 }
0030
0031 void CreateSurveyRcds::analyze(const edm::Event& event, const edm::EventSetup& setup) {
0032
0033 const TrackerTopology* const tTopo = &setup.getData(tTopoToken_);
0034 const GeometricDet* geom = &setup.getData(geomDetToken_);
0035 const PTrackerParameters& ptp = setup.getData(ptpToken_);
0036 const PTrackerAdditionalParametersPerDet* ptitp = &setup.getData(ptitpToken_);
0037 TrackerGeometry* tracker = TrackerGeomBuilderFromGeometricDet().build(geom, ptitp, ptp, tTopo);
0038
0039
0040 if (m_inputGeom == "sqlite") {
0041
0042 const Alignments* alignments = &setup.getData(aliToken_);
0043 const AlignmentErrorsExtended* alignmentErrors = &setup.getData(aliErrToken_);
0044
0045
0046 GeometryAligner aligner;
0047 aligner.applyAlignments<TrackerGeometry>(&(*tracker), alignments, alignmentErrors, AlignTransform());
0048 }
0049
0050 addComponent(new AlignableTracker(tracker, tTopo));
0051
0052 Alignable* ali = detector();
0053 if (m_inputGeom == "generated") {
0054 setGeometry(ali);
0055 }
0056
0057 setSurveyErrors(ali);
0058 }
0059
0060 void CreateSurveyRcds::setGeometry(Alignable* ali) {
0061 const align::Alignables& comp = ali->components();
0062 unsigned int nComp = comp.size();
0063
0064
0065 bool usecomps = true;
0066 if ((ali->alignableObjectId() == 2) && (nComp >= 1))
0067 usecomps = false;
0068 for (unsigned int i = 0; i < nComp; ++i) {
0069 if (usecomps)
0070 setGeometry(comp[i]);
0071 }
0072 DetId id(ali->id());
0073 int subdetlevel = id.subdetId();
0074 int level = ali->alignableObjectId();
0075
0076
0077 if (m_generatedRandom) {
0078 if (subdetlevel > 0) {
0079 AlgebraicVector value = getStructureErrors(level, subdetlevel);
0080
0081 double value0 = CLHEP::RandGauss::shoot(0, value[0]);
0082 double value1 = CLHEP::RandGauss::shoot(0, value[1]);
0083 double value2 = CLHEP::RandGauss::shoot(0, value[2]);
0084 double value3 = CLHEP::RandGauss::shoot(0, value[3]);
0085 double value4 = CLHEP::RandGauss::shoot(0, value[4]);
0086 double value5 = CLHEP::RandGauss::shoot(0, value[5]);
0087
0088
0089 align::LocalVector diffR(value0, value1, value2);
0090 align::Scalar diffWx = value3;
0091 align::Scalar diffWy = value4;
0092 align::Scalar diffWz = value5;
0093 ali->move(ali->surface().toGlobal(diffR));
0094 ali->rotateAroundLocalX(diffWx);
0095 ali->rotateAroundLocalY(diffWy);
0096 ali->rotateAroundLocalZ(diffWz);
0097 }
0098 }
0099
0100
0101 if (m_generatedSimple) {
0102 if ((level == 2) || ((level == 1) && (ali->mother()->alignableObjectId() != 2))) {
0103 const double constMis = m_inputSimpleMis;
0104 const double dAngle = constMis / ali->surface().length();
0105
0106 double value0 = CLHEP::RandGauss::shoot(0, constMis);
0107 double value1 = CLHEP::RandGauss::shoot(0, constMis);
0108 double value2 = CLHEP::RandGauss::shoot(0, constMis);
0109 double value3 = CLHEP::RandGauss::shoot(0, dAngle);
0110 double value4 = CLHEP::RandGauss::shoot(0, dAngle);
0111 double value5 = CLHEP::RandGauss::shoot(0, dAngle);
0112
0113 align::LocalVector diffR(value0, value1, value2);
0114 ali->move(ali->surface().toGlobal(diffR));
0115 align::Scalar diffWx = value3;
0116 align::Scalar diffWy = value4;
0117 align::Scalar diffWz = value5;
0118 ali->rotateAroundLocalX(diffWx);
0119 ali->rotateAroundLocalY(diffWy);
0120 ali->rotateAroundLocalZ(diffWz);
0121 }
0122 }
0123 }
0124
0125 void CreateSurveyRcds::setSurveyErrors(Alignable* ali) {
0126 const align::Alignables& comp = ali->components();
0127 unsigned int nComp = comp.size();
0128
0129
0130 for (unsigned int i = 0; i < nComp; ++i) {
0131 setSurveyErrors(comp[i]);
0132 }
0133
0134 DetId id(ali->id());
0135 int subdetlevel = id.subdetId();
0136 int level = ali->alignableObjectId();
0137
0138 AlgebraicVector error = getStructureErrors(level, subdetlevel);
0139
0140 double error0 = error[0];
0141 double error1 = error[1];
0142 double error2 = error[2];
0143 double error3 = error[3];
0144 double error4 = error[4];
0145 double error5 = error[5];
0146
0147
0148
0149
0150 if ((level <= 2) && (subdetlevel == 3)) {
0151 error0 = 0.01;
0152 error5 = 0.001;
0153 if ((level == 2) && (nComp == 2)) {
0154 error1 = 0.01;
0155 }
0156 }
0157
0158 if ((level <= 2) && (subdetlevel == 4)) {
0159 error0 = 0.01;
0160 error1 = 0.01;
0161 error2 = 0.01;
0162 error3 = 0.001;
0163 error4 = 0.001;
0164 error5 = 0.001;
0165
0166
0167
0168
0169 }
0170 if ((level == 23) && (subdetlevel == 4)) {
0171 error0 = 0.02;
0172 error1 = 0.02;
0173 error2 = 0.03;
0174 error3 = 0.0002;
0175 error4 = 0.0002;
0176 error5 = 0.0002;
0177 }
0178 if ((level == 22) && (subdetlevel == 4)) {
0179 error0 = 0.01;
0180 error1 = 0.01;
0181 error2 = 0.01;
0182 error3 = 0.0002;
0183 error4 = 0.0002;
0184 error5 = 0.0002;
0185 }
0186
0187 if ((level <= 2) && (subdetlevel == 5)) {
0188
0189 error0 = 0.015;
0190 error1 = 0.003;
0191 error2 = 0.003;
0192 error3 = 0.0002;
0193 error4 = 0.0002;
0194 error5 = 0.001;
0195 if ((level == 2) && (nComp == 2)) {
0196 error1 = 0.015;
0197 }
0198 }
0199 if ((level == 27) && (subdetlevel == 5)) {
0200 error0 = 0.02;
0201 error1 = 0.02;
0202 error2 = 0.03;
0203 error3 = 0.001;
0204 error4 = 0.001;
0205 error5 = 0.001;
0206 }
0207
0208 if ((level <= 2) && (subdetlevel == 6)) {
0209 error0 = 0.02;
0210 error5 = 0.0005;
0211 if ((level == 2) && (nComp == 2)) {
0212 error1 = 0.02;
0213 }
0214 }
0215 if ((level == 34) && (subdetlevel == 6)) {
0216 error0 = 0.01;
0217 error1 = 0.01;
0218 error2 = 0.02;
0219 error3 = 0.00005;
0220 error4 = 0.00005;
0221 error5 = 0.00005;
0222 }
0223 if ((level == 33) && (subdetlevel == 6)) {
0224 error0 = 0.01;
0225 error1 = 0.01;
0226 error2 = 0.02;
0227 error3 = 0.0001;
0228 error4 = 0.0001;
0229 error5 = 0.0001;
0230 }
0231 if ((level == 32) && (subdetlevel == 6)) {
0232 error0 = 0.007;
0233 error1 = 0.007;
0234 error2 = 0.015;
0235 error3 = 0.00015;
0236 error4 = 0.00015;
0237 error5 = 0.00015;
0238 }
0239
0240
0241
0242 align::ErrorMatrix error_Matrix;
0243 double* errorData = error_Matrix.Array();
0244 errorData[0] = error0 * error0;
0245 errorData[2] = error1 * error1;
0246 errorData[5] = error2 * error2;
0247 errorData[9] = error3 * error3;
0248 errorData[14] = error4 * error4;
0249 errorData[20] = error5 * error5;
0250 errorData[1] = 0.0;
0251 errorData[3] = 0.0;
0252 errorData[4] = 0.0;
0253 errorData[6] = 0.0;
0254 errorData[7] = 0.0;
0255 errorData[8] = 0.0;
0256 errorData[10] = 0.0;
0257 errorData[11] = 0.0;
0258 errorData[12] = 0.0;
0259 errorData[13] = 0.0;
0260 errorData[15] = 0.0;
0261 errorData[16] = 0.0;
0262 errorData[17] = 0.0;
0263 errorData[18] = 0.0;
0264 errorData[19] = 0.0;
0265
0266 ali->setSurvey(new SurveyDet(ali->surface(), error_Matrix));
0267 }
0268
0269
0270
0271
0272 AlgebraicVector CreateSurveyRcds::getStructurePlacements(int level, int subdetlevel) {
0273 AlgebraicVector deltaRW(6);
0274 deltaRW(1) = 0.0;
0275 deltaRW(2) = 0.0;
0276 deltaRW(3) = 0.0;
0277 deltaRW(4) = 0.0;
0278 deltaRW(5) = 0.0;
0279 deltaRW(6) = 0.0;
0280
0281 if ((level == 37) && (subdetlevel == 1)) {
0282 deltaRW(1) = 0.3;
0283 deltaRW(2) = 0.3;
0284 deltaRW(3) = 0.3;
0285 deltaRW(4) = 0.0017;
0286 deltaRW(5) = 0.0017;
0287 deltaRW(6) = 0.0017;
0288 }
0289
0290 if ((level == 38) && (subdetlevel == 3)) {
0291 deltaRW(1) = 0.3;
0292 deltaRW(2) = 0.3;
0293 deltaRW(3) = 0.3;
0294 deltaRW(4) = 0.0004;
0295 deltaRW(5) = 0.0004;
0296 deltaRW(6) = 0.0004;
0297 }
0298
0299 if ((level == 39) && (subdetlevel == 1)) {
0300 deltaRW(1) = 0.0;
0301 deltaRW(2) = 0.0;
0302 deltaRW(3) = 0.0;
0303 deltaRW(4) = 0.0;
0304 deltaRW(5) = 0.0;
0305 deltaRW(6) = 0.0;
0306 }
0307
0308 if ((level == 7) && (subdetlevel == 1)) {
0309 deltaRW(1) = 0.2;
0310 deltaRW(2) = 0.2;
0311 deltaRW(3) = 0.2;
0312 deltaRW(4) = 0.003;
0313 deltaRW(5) = 0.003;
0314 deltaRW(6) = 0.003;
0315 }
0316 if ((level == 6) && (subdetlevel == 1)) {
0317 deltaRW(1) = 0.05;
0318 deltaRW(2) = 0.05;
0319 deltaRW(3) = 0.05;
0320 deltaRW(4) = 0.0008;
0321 deltaRW(5) = 0.0008;
0322 deltaRW(6) = 0.0008;
0323 }
0324 if ((level == 5) && (subdetlevel == 1)) {
0325 deltaRW(1) = 0.02;
0326 deltaRW(2) = 0.02;
0327 deltaRW(3) = 0.02;
0328 deltaRW(4) = 0.0004;
0329 deltaRW(5) = 0.0004;
0330 deltaRW(6) = 0.0004;
0331 }
0332 if ((level == 4) && (subdetlevel == 1)) {
0333 deltaRW(1) = 0.01;
0334 deltaRW(2) = 0.01;
0335 deltaRW(3) = 0.005;
0336 deltaRW(4) = 0.0002;
0337 deltaRW(5) = 0.0002;
0338 deltaRW(6) = 0.0002;
0339 }
0340 if ((level == 2) && (subdetlevel == 1)) {
0341 deltaRW(1) = 0.005;
0342 deltaRW(2) = 0.005;
0343 deltaRW(3) = 0.003;
0344 deltaRW(4) = 0.001;
0345 deltaRW(5) = 0.001;
0346 deltaRW(6) = 0.001;
0347 }
0348 if ((level == 1) && (subdetlevel == 1)) {
0349 deltaRW(1) = 0.005;
0350 deltaRW(2) = 0.005;
0351 deltaRW(3) = 0.003;
0352 deltaRW(4) = 0.001;
0353 deltaRW(5) = 0.001;
0354 deltaRW(6) = 0.001;
0355 }
0356
0357 if ((level == 13) && (subdetlevel == 2)) {
0358 deltaRW(1) = 0.2;
0359 deltaRW(2) = 0.2;
0360 deltaRW(3) = 0.2;
0361 deltaRW(4) = 0.0017;
0362 deltaRW(5) = 0.0017;
0363 deltaRW(6) = 0.0017;
0364 }
0365 if ((level == 12) && (subdetlevel == 2)) {
0366 deltaRW(1) = 0.05;
0367 deltaRW(2) = 0.05;
0368 deltaRW(3) = 0.05;
0369 deltaRW(4) = 0.0004;
0370 deltaRW(5) = 0.0004;
0371 deltaRW(6) = 0.0004;
0372 }
0373 if ((level == 11) && (subdetlevel == 2)) {
0374 deltaRW(1) = 0.02;
0375 deltaRW(2) = 0.02;
0376 deltaRW(3) = 0.02;
0377 deltaRW(4) = 0.001;
0378 deltaRW(5) = 0.001;
0379 deltaRW(6) = 0.001;
0380 }
0381 if ((level == 10) && (subdetlevel == 2)) {
0382 deltaRW(1) = 0.01;
0383 deltaRW(2) = 0.01;
0384 deltaRW(3) = 0.01;
0385 deltaRW(4) = 0.001;
0386 deltaRW(5) = 0.001;
0387 deltaRW(6) = 0.001;
0388 }
0389 if ((level == 9) && (subdetlevel == 2)) {
0390 deltaRW(1) = 0.01;
0391 deltaRW(2) = 0.01;
0392 deltaRW(3) = 0.005;
0393 deltaRW(4) = 0.002;
0394 deltaRW(5) = 0.002;
0395 deltaRW(6) = 0.002;
0396 }
0397 if ((level == 2) && (subdetlevel == 2)) {
0398 deltaRW(1) = 0.005;
0399 deltaRW(2) = 0.005;
0400 deltaRW(3) = 0.003;
0401 deltaRW(4) = 0.001;
0402 deltaRW(5) = 0.001;
0403 deltaRW(6) = 0.001;
0404 }
0405 if ((level == 1) && (subdetlevel == 2)) {
0406 deltaRW(1) = 0.005;
0407 deltaRW(2) = 0.005;
0408 deltaRW(3) = 0.003;
0409 deltaRW(4) = 0.001;
0410 deltaRW(5) = 0.001;
0411 deltaRW(6) = 0.001;
0412 }
0413
0414 if ((level == 20) && (subdetlevel == 3)) {
0415 deltaRW(1) = 0.2;
0416 deltaRW(2) = 0.2;
0417 deltaRW(3) = 0.2;
0418 deltaRW(4) = 0.0017;
0419 deltaRW(5) = 0.0017;
0420 deltaRW(6) = 0.0017;
0421 }
0422 if ((level == 19) && (subdetlevel == 3)) {
0423 deltaRW(1) = 0.1;
0424 deltaRW(2) = 0.1;
0425 deltaRW(3) = 0.1;
0426 deltaRW(4) = 0.0008;
0427 deltaRW(5) = 0.0008;
0428 deltaRW(6) = 0.0008;
0429 }
0430 if ((level == 18) && (subdetlevel == 3)) {
0431 deltaRW(1) = 0.04;
0432 deltaRW(2) = 0.04;
0433 deltaRW(3) = 0.02;
0434 deltaRW(4) = 0.0006;
0435 deltaRW(5) = 0.0006;
0436 deltaRW(6) = 0.0006;
0437 }
0438 if ((level == 17) && (subdetlevel == 3)) {
0439 deltaRW(1) = 0.03;
0440 deltaRW(2) = 0.03;
0441 deltaRW(3) = 0.015;
0442 deltaRW(4) = 0.0004;
0443 deltaRW(5) = 0.0004;
0444 deltaRW(6) = 0.0004;
0445 }
0446 if ((level == 16) && (subdetlevel == 3)) {
0447 deltaRW(1) = 0.01;
0448 deltaRW(2) = 0.01;
0449 deltaRW(3) = 0.01;
0450 deltaRW(4) = 0.0004;
0451 deltaRW(5) = 0.0002;
0452 deltaRW(6) = 0.0002;
0453 }
0454 if ((level == 15) && (subdetlevel == 3)) {
0455 deltaRW(1) = 0.01;
0456 deltaRW(2) = 0.01;
0457 deltaRW(3) = 0.01;
0458 deltaRW(4) = 0.0004;
0459 deltaRW(5) = 0.0002;
0460 deltaRW(6) = 0.0002;
0461 }
0462 if ((level == 2) && (subdetlevel == 3)) {
0463 deltaRW(1) = 0.005;
0464 deltaRW(2) = 0.005;
0465 deltaRW(3) = 0.005;
0466 deltaRW(4) = 0.001;
0467 deltaRW(5) = 0.0005;
0468 deltaRW(6) = 0.0005;
0469 }
0470 if ((level == 1) && (subdetlevel == 3)) {
0471 deltaRW(1) = 0.005;
0472 deltaRW(2) = 0.005;
0473 deltaRW(3) = 0.005;
0474 deltaRW(4) = 0.001;
0475 deltaRW(5) = 0.0005;
0476 deltaRW(6) = 0.0005;
0477 }
0478
0479 if ((level == 25) && (subdetlevel == 4)) {
0480 deltaRW(1) = 0.2;
0481 deltaRW(2) = 0.2;
0482 deltaRW(3) = 0.2;
0483 deltaRW(4) = 0.0013;
0484 deltaRW(5) = 0.0013;
0485 deltaRW(6) = 0.0013;
0486 }
0487 if ((level == 24) && (subdetlevel == 4)) {
0488 deltaRW(1) = 0.05;
0489 deltaRW(2) = 0.05;
0490 deltaRW(3) = 0.05;
0491 deltaRW(4) = 0.0004;
0492 deltaRW(5) = 0.0004;
0493 deltaRW(6) = 0.0004;
0494 }
0495 if ((level == 23) && (subdetlevel == 4)) {
0496 deltaRW(1) = 0.01;
0497 deltaRW(2) = 0.01;
0498 deltaRW(3) = 0.01;
0499 deltaRW(4) = 0.0001;
0500 deltaRW(5) = 0.0001;
0501 deltaRW(6) = 0.0001;
0502 }
0503 if ((level == 22) && (subdetlevel == 4)) {
0504 deltaRW(1) = 0.005;
0505 deltaRW(2) = 0.005;
0506 deltaRW(3) = 0.005;
0507 deltaRW(4) = 0.0001;
0508 deltaRW(5) = 0.0001;
0509 deltaRW(6) = 0.0001;
0510 }
0511 if ((level == 2) && (subdetlevel == 4)) {
0512 deltaRW(1) = 0.005;
0513 deltaRW(2) = 0.005;
0514 deltaRW(3) = 0.005;
0515 deltaRW(4) = 0.0005;
0516 deltaRW(5) = 0.0005;
0517 deltaRW(6) = 0.0005;
0518 }
0519 if ((level == 1) && (subdetlevel == 4)) {
0520 deltaRW(1) = 0.005;
0521 deltaRW(2) = 0.005;
0522 deltaRW(3) = 0.005;
0523 deltaRW(4) = 0.0005;
0524 deltaRW(5) = 0.0005;
0525 deltaRW(6) = 0.0005;
0526 }
0527
0528 if ((level == 30) && (subdetlevel == 5)) {
0529 deltaRW(1) = 0.2;
0530 deltaRW(2) = 0.2;
0531 deltaRW(3) = 0.2;
0532 deltaRW(4) = 0.0008;
0533 deltaRW(5) = 0.0008;
0534 deltaRW(6) = 0.0008;
0535 }
0536 if ((level == 29) && (subdetlevel == 5)) {
0537 deltaRW(1) = 0.014;
0538 deltaRW(2) = 0.014;
0539 deltaRW(3) = 0.05;
0540 deltaRW(4) = 0.0001;
0541 deltaRW(5) = 0.0001;
0542 deltaRW(6) = 0.0001;
0543 }
0544 if ((level == 28) && (subdetlevel == 5)) {
0545 deltaRW(1) = 0.02;
0546 deltaRW(2) = 0.02;
0547 deltaRW(3) = 0.02;
0548 deltaRW(4) = 0.0001;
0549 deltaRW(5) = 0.0001;
0550 deltaRW(6) = 0.0001;
0551 }
0552 if ((level == 27) && (subdetlevel == 5)) {
0553 deltaRW(1) = 0.01;
0554 deltaRW(2) = 0.01;
0555 deltaRW(3) = 0.02;
0556 deltaRW(4) = 0.0001;
0557 deltaRW(5) = 0.0001;
0558 deltaRW(6) = 0.0001;
0559 }
0560 if ((level == 2) && (subdetlevel == 5)) {
0561 deltaRW(1) = 0.003;
0562 deltaRW(2) = 0.003;
0563 deltaRW(3) = 0.01;
0564 deltaRW(4) = 0.0002;
0565 deltaRW(5) = 0.0002;
0566 deltaRW(6) = 0.0002;
0567 }
0568 if ((level == 1) && (subdetlevel == 5)) {
0569 deltaRW(1) = 0.003;
0570 deltaRW(2) = 0.003;
0571 deltaRW(3) = 0.01;
0572 deltaRW(4) = 0.0002;
0573 deltaRW(5) = 0.0002;
0574 deltaRW(6) = 0.0002;
0575 }
0576
0577 if ((level == 36) && (subdetlevel == 6)) {
0578 deltaRW(1) = 0.2;
0579 deltaRW(2) = 0.2;
0580 deltaRW(3) = 0.2;
0581 deltaRW(4) = 0.0008;
0582 deltaRW(5) = 0.0008;
0583 deltaRW(6) = 0.0008;
0584 }
0585 if ((level == 35) && (subdetlevel == 6)) {
0586 deltaRW(1) = 0.05;
0587 deltaRW(2) = 0.05;
0588 deltaRW(3) = 0.05;
0589 deltaRW(4) = 0.0003;
0590 deltaRW(5) = 0.0003;
0591 deltaRW(6) = 0.0003;
0592 }
0593 if ((level == 34) && (subdetlevel == 6)) {
0594 deltaRW(1) = 0.01;
0595 deltaRW(2) = 0.01;
0596 deltaRW(3) = 0.02;
0597 deltaRW(4) = 0.00005;
0598 deltaRW(5) = 0.00005;
0599 deltaRW(6) = 0.00005;
0600 }
0601 if ((level == 33) && (subdetlevel == 6)) {
0602 deltaRW(1) = 0.01;
0603 deltaRW(2) = 0.01;
0604 deltaRW(3) = 0.02;
0605 deltaRW(4) = 0.0001;
0606 deltaRW(5) = 0.0001;
0607 deltaRW(6) = 0.0001;
0608 }
0609 if ((level == 32) && (subdetlevel == 6)) {
0610 deltaRW(1) = 0.007;
0611 deltaRW(2) = 0.007;
0612 deltaRW(3) = 0.015;
0613 deltaRW(4) = 0.00015;
0614 deltaRW(5) = 0.00015;
0615 deltaRW(6) = 0.00015;
0616 }
0617 if ((level == 2) && (subdetlevel == 6)) {
0618 deltaRW(1) = 0.002;
0619 deltaRW(2) = 0.002;
0620 deltaRW(3) = 0.005;
0621 deltaRW(4) = 0.0001;
0622 deltaRW(5) = 0.0001;
0623 deltaRW(6) = 0.0001;
0624 }
0625 if ((level == 1) && (subdetlevel == 6)) {
0626 deltaRW(1) = 0.002;
0627 deltaRW(2) = 0.002;
0628 deltaRW(3) = 0.005;
0629 deltaRW(4) = 0.0001;
0630 deltaRW(5) = 0.0001;
0631 deltaRW(6) = 0.0001;
0632 }
0633
0634 return deltaRW;
0635 }
0636
0637
0638
0639 AlgebraicVector CreateSurveyRcds::getStructureErrors(int level, int subdetlevel) {
0640 AlgebraicVector deltaRW(6);
0641 deltaRW(1) = 0.0;
0642 deltaRW(2) = 0.0;
0643 deltaRW(3) = 0.0;
0644 deltaRW(4) = 0.0;
0645 deltaRW(5) = 0.0;
0646 deltaRW(6) = 0.0;
0647
0648 if ((level == 37) && (subdetlevel == 1)) {
0649 deltaRW(1) = 0.2;
0650 deltaRW(2) = 0.2;
0651 deltaRW(3) = 0.2;
0652 deltaRW(4) = 0.0017;
0653 deltaRW(5) = 0.0017;
0654 deltaRW(6) = 0.0017;
0655 }
0656
0657 if ((level == 38) && (subdetlevel == 3)) {
0658 deltaRW(1) = 0.2;
0659 deltaRW(2) = 0.2;
0660 deltaRW(3) = 0.2;
0661 deltaRW(4) = 0.0004;
0662 deltaRW(5) = 0.0004;
0663 deltaRW(6) = 0.0004;
0664 }
0665
0666 if ((level == 39) && (subdetlevel == 1)) {
0667 deltaRW(1) = 0.0;
0668 deltaRW(2) = 0.0;
0669 deltaRW(3) = 0.0;
0670 deltaRW(4) = 0.0;
0671 deltaRW(5) = 0.0;
0672 deltaRW(6) = 0.0;
0673 }
0674
0675 if ((level == 7) && (subdetlevel == 1)) {
0676 deltaRW(1) = 0.05;
0677 deltaRW(2) = 0.05;
0678 deltaRW(3) = 0.1;
0679 deltaRW(4) = 0.0008;
0680 deltaRW(5) = 0.0008;
0681 deltaRW(6) = 0.0008;
0682 }
0683 if ((level == 6) && (subdetlevel == 1)) {
0684 deltaRW(1) = 0.015;
0685 deltaRW(2) = 0.015;
0686 deltaRW(3) = 0.03;
0687 deltaRW(4) = 0.0003;
0688 deltaRW(5) = 0.0003;
0689 deltaRW(6) = 0.0003;
0690 }
0691 if ((level == 5) && (subdetlevel == 1)) {
0692 deltaRW(1) = 0.005;
0693 deltaRW(2) = 0.005;
0694 deltaRW(3) = 0.01;
0695 deltaRW(4) = 0.0001;
0696 deltaRW(5) = 0.0001;
0697 deltaRW(6) = 0.0001;
0698 }
0699 if ((level == 4) && (subdetlevel == 1)) {
0700 deltaRW(1) = 0.001;
0701 deltaRW(2) = 0.001;
0702 deltaRW(3) = 0.002;
0703 deltaRW(4) = 0.00005;
0704 deltaRW(5) = 0.00005;
0705 deltaRW(6) = 0.00005;
0706 }
0707 if ((level == 2) && (subdetlevel == 1)) {
0708 deltaRW(1) = 0.0005;
0709 deltaRW(2) = 0.001;
0710 deltaRW(3) = 0.001;
0711 deltaRW(4) = 0.0001;
0712 deltaRW(5) = 0.0001;
0713 deltaRW(6) = 0.0001;
0714 }
0715 if ((level == 1) && (subdetlevel == 1)) {
0716 deltaRW(1) = 0.0005;
0717 deltaRW(2) = 0.001;
0718 deltaRW(3) = 0.001;
0719 deltaRW(4) = 0.0001;
0720 deltaRW(5) = 0.0001;
0721 deltaRW(6) = 0.0001;
0722 }
0723
0724 if ((level == 13) && (subdetlevel == 2)) {
0725 deltaRW(1) = 0.05;
0726 deltaRW(2) = 0.05;
0727 deltaRW(3) = 0.1;
0728 deltaRW(4) = 0.0004;
0729 deltaRW(5) = 0.0004;
0730 deltaRW(6) = 0.0004;
0731 }
0732 if ((level == 12) && (subdetlevel == 2)) {
0733 deltaRW(1) = 0.015;
0734 deltaRW(2) = 0.015;
0735 deltaRW(3) = 0.03;
0736 deltaRW(4) = 0.00015;
0737 deltaRW(5) = 0.00015;
0738 deltaRW(6) = 0.00015;
0739 }
0740 if ((level == 11) && (subdetlevel == 2)) {
0741 deltaRW(1) = 0.005;
0742 deltaRW(2) = 0.005;
0743 deltaRW(3) = 0.01;
0744 deltaRW(4) = 0.0001;
0745 deltaRW(5) = 0.0001;
0746 deltaRW(6) = 0.0001;
0747 }
0748 if ((level == 10) && (subdetlevel == 2)) {
0749 deltaRW(1) = 0.001;
0750 deltaRW(2) = 0.001;
0751 deltaRW(3) = 0.002;
0752 deltaRW(4) = 0.0001;
0753 deltaRW(5) = 0.0001;
0754 deltaRW(6) = 0.0001;
0755 }
0756 if ((level == 9) && (subdetlevel == 2)) {
0757 deltaRW(1) = 0.001;
0758 deltaRW(2) = 0.0008;
0759 deltaRW(3) = 0.0006;
0760 deltaRW(4) = 0.0002;
0761 deltaRW(5) = 0.0002;
0762 deltaRW(6) = 0.0002;
0763 }
0764 if ((level == 2) && (subdetlevel == 2)) {
0765 deltaRW(1) = 0.0005;
0766 deltaRW(2) = 0.0004;
0767 deltaRW(3) = 0.0006;
0768 deltaRW(4) = 0.0001;
0769 deltaRW(5) = 0.0003;
0770 deltaRW(6) = 0.0001;
0771 }
0772 if ((level == 1) && (subdetlevel == 2)) {
0773 deltaRW(1) = 0.0005;
0774 deltaRW(2) = 0.0004;
0775 deltaRW(3) = 0.0006;
0776 deltaRW(4) = 0.0001;
0777 deltaRW(5) = 0.0003;
0778 deltaRW(6) = 0.0001;
0779 }
0780
0781 if ((level == 20) && (subdetlevel == 3)) {
0782 deltaRW(1) = 0.08;
0783 deltaRW(2) = 0.08;
0784 deltaRW(3) = 0.04;
0785 deltaRW(4) = 0.0017;
0786 deltaRW(5) = 0.0017;
0787 deltaRW(6) = 0.0017;
0788 }
0789 if ((level == 19) && (subdetlevel == 3)) {
0790 deltaRW(1) = 0.04;
0791 deltaRW(2) = 0.04;
0792 deltaRW(3) = 0.02;
0793 deltaRW(4) = 0.0003;
0794 deltaRW(5) = 0.0003;
0795 deltaRW(6) = 0.0003;
0796 }
0797 if ((level == 18) && (subdetlevel == 3)) {
0798 deltaRW(1) = 0.02;
0799 deltaRW(2) = 0.02;
0800 deltaRW(3) = 0.01;
0801 deltaRW(4) = 0.0006;
0802 deltaRW(5) = 0.0006;
0803 deltaRW(6) = 0.0006;
0804 }
0805 if ((level == 17) && (subdetlevel == 3)) {
0806 deltaRW(1) = 0.01;
0807 deltaRW(2) = 0.01;
0808 deltaRW(3) = 0.005;
0809 deltaRW(4) = 0.0002;
0810 deltaRW(5) = 0.0002;
0811 deltaRW(6) = 0.0002;
0812 }
0813 if ((level == 16) && (subdetlevel == 3)) {
0814 deltaRW(1) = 0.004;
0815 deltaRW(2) = 0.004;
0816 deltaRW(3) = 0.008;
0817 deltaRW(4) = 0.0002;
0818 deltaRW(5) = 0.0001;
0819 deltaRW(6) = 0.0001;
0820 }
0821 if ((level == 15) && (subdetlevel == 3)) {
0822 deltaRW(1) = 0.004;
0823 deltaRW(2) = 0.004;
0824 deltaRW(3) = 0.008;
0825 deltaRW(4) = 0.0002;
0826 deltaRW(5) = 0.0001;
0827 deltaRW(6) = 0.0001;
0828 }
0829 if ((level == 2) && (subdetlevel == 3)) {
0830 deltaRW(1) = 0.002;
0831 deltaRW(2) = 0.002;
0832 deltaRW(3) = 0.004;
0833 deltaRW(4) = 0.0004;
0834 deltaRW(5) = 0.0002;
0835 deltaRW(6) = 0.0002;
0836 }
0837 if ((level == 1) && (subdetlevel == 3)) {
0838 deltaRW(1) = 0.002;
0839 deltaRW(2) = 0.002;
0840 deltaRW(3) = 0.004;
0841 deltaRW(4) = 0.0004;
0842 deltaRW(5) = 0.0002;
0843 deltaRW(6) = 0.0002;
0844 }
0845
0846 if ((level == 25) && (subdetlevel == 4)) {
0847 deltaRW(1) = 0.05;
0848 deltaRW(2) = 0.05;
0849 deltaRW(3) = 0.1;
0850 deltaRW(4) = 0.0003;
0851 deltaRW(5) = 0.0003;
0852 deltaRW(6) = 0.0003;
0853 }
0854 if ((level == 24) && (subdetlevel == 4)) {
0855 deltaRW(1) = 0.01;
0856 deltaRW(2) = 0.01;
0857 deltaRW(3) = 0.02;
0858 deltaRW(4) = 0.0001;
0859 deltaRW(5) = 0.0001;
0860 deltaRW(6) = 0.0001;
0861 }
0862 if ((level == 23) && (subdetlevel == 4)) {
0863 deltaRW(1) = 0.004;
0864 deltaRW(2) = 0.004;
0865 deltaRW(3) = 0.005;
0866 deltaRW(4) = 0.00004;
0867 deltaRW(5) = 0.00004;
0868 deltaRW(6) = 0.00004;
0869 }
0870 if ((level == 22) && (subdetlevel == 4)) {
0871 deltaRW(1) = 0.002;
0872 deltaRW(2) = 0.002;
0873 deltaRW(3) = 0.002;
0874 deltaRW(4) = 0.00004;
0875 deltaRW(5) = 0.00004;
0876 deltaRW(6) = 0.00004;
0877 }
0878 if ((level == 2) && (subdetlevel == 4)) {
0879 deltaRW(1) = 0.002;
0880 deltaRW(2) = 0.002;
0881 deltaRW(3) = 0.002;
0882 deltaRW(4) = 0.0002;
0883 deltaRW(5) = 0.0002;
0884 deltaRW(6) = 0.0002;
0885 }
0886 if ((level == 1) && (subdetlevel == 4)) {
0887 deltaRW(1) = 0.002;
0888 deltaRW(2) = 0.002;
0889 deltaRW(3) = 0.002;
0890 deltaRW(4) = 0.0002;
0891 deltaRW(5) = 0.0002;
0892 deltaRW(6) = 0.0002;
0893 }
0894
0895 if ((level == 30) && (subdetlevel == 5)) {
0896 deltaRW(1) = 0.06;
0897 deltaRW(2) = 0.06;
0898 deltaRW(3) = 0.06;
0899 deltaRW(4) = 0.00025;
0900 deltaRW(5) = 0.00025;
0901 deltaRW(6) = 0.00025;
0902 }
0903 if ((level == 29) && (subdetlevel == 5)) {
0904 deltaRW(1) = 0.014;
0905 deltaRW(2) = 0.014;
0906 deltaRW(3) = 0.05;
0907 deltaRW(4) = 0.0001;
0908 deltaRW(5) = 0.0001;
0909 deltaRW(6) = 0.0001;
0910 }
0911 if ((level == 28) && (subdetlevel == 5)) {
0912 deltaRW(1) = 0.02;
0913 deltaRW(2) = 0.02;
0914 deltaRW(3) = 0.02;
0915 deltaRW(4) = 0.0001;
0916 deltaRW(5) = 0.0001;
0917 deltaRW(6) = 0.0001;
0918 }
0919 if ((level == 27) && (subdetlevel == 5)) {
0920 deltaRW(1) = 0.01;
0921 deltaRW(2) = 0.01;
0922 deltaRW(3) = 0.02;
0923 deltaRW(4) = 0.0001;
0924 deltaRW(5) = 0.0001;
0925 deltaRW(6) = 0.0001;
0926 }
0927 if ((level == 2) && (subdetlevel == 5)) {
0928 deltaRW(1) = 0.003;
0929 deltaRW(2) = 0.003;
0930 deltaRW(3) = 0.01;
0931 deltaRW(4) = 0.0002;
0932 deltaRW(5) = 0.0002;
0933 deltaRW(6) = 0.0002;
0934 }
0935 if ((level == 1) && (subdetlevel == 5)) {
0936 deltaRW(1) = 0.003;
0937 deltaRW(2) = 0.003;
0938 deltaRW(3) = 0.01;
0939 deltaRW(4) = 0.0002;
0940 deltaRW(5) = 0.0002;
0941 deltaRW(6) = 0.0002;
0942 }
0943
0944 if ((level == 36) && (subdetlevel == 6)) {
0945 deltaRW(1) = 0.06;
0946 deltaRW(2) = 0.06;
0947 deltaRW(3) = 0.1;
0948 deltaRW(4) = 0.0003;
0949 deltaRW(5) = 0.0003;
0950 deltaRW(6) = 0.0003;
0951 }
0952 if ((level == 35) && (subdetlevel == 6)) {
0953 deltaRW(1) = 0.015;
0954 deltaRW(2) = 0.015;
0955 deltaRW(3) = 0.03;
0956 deltaRW(4) = 0.0001;
0957 deltaRW(5) = 0.0001;
0958 deltaRW(6) = 0.0001;
0959 }
0960 if ((level == 34) && (subdetlevel == 6)) {
0961 deltaRW(1) = 0.01;
0962 deltaRW(2) = 0.01;
0963 deltaRW(3) = 0.02;
0964 deltaRW(4) = 0.00005;
0965 deltaRW(5) = 0.00005;
0966 deltaRW(6) = 0.00005;
0967 }
0968 if ((level == 33) && (subdetlevel == 6)) {
0969 deltaRW(1) = 0.01;
0970 deltaRW(2) = 0.01;
0971 deltaRW(3) = 0.02;
0972 deltaRW(4) = 0.0001;
0973 deltaRW(5) = 0.0001;
0974 deltaRW(6) = 0.0001;
0975 }
0976 if ((level == 32) && (subdetlevel == 6)) {
0977 deltaRW(1) = 0.007;
0978 deltaRW(2) = 0.007;
0979 deltaRW(3) = 0.015;
0980 deltaRW(4) = 0.00015;
0981 deltaRW(5) = 0.00015;
0982 deltaRW(6) = 0.00015;
0983 }
0984 if ((level == 2) && (subdetlevel == 6)) {
0985 deltaRW(1) = 0.002;
0986 deltaRW(2) = 0.002;
0987 deltaRW(3) = 0.005;
0988 deltaRW(4) = 0.0001;
0989 deltaRW(5) = 0.0001;
0990 deltaRW(6) = 0.0001;
0991 }
0992 if ((level == 1) && (subdetlevel == 6)) {
0993 deltaRW(1) = 0.002;
0994 deltaRW(2) = 0.002;
0995 deltaRW(3) = 0.005;
0996 deltaRW(4) = 0.0001;
0997 deltaRW(5) = 0.0001;
0998 deltaRW(6) = 0.0001;
0999 }
1000
1001 return deltaRW;
1002 }
1003
1004
1005
1006 #include "FWCore/Framework/interface/MakerMacros.h"
1007
1008 DEFINE_FWK_MODULE(CreateSurveyRcds);