68 #include "TStopwatch.h" 204 m_muonCollectionTag(
cfg.getParameter<
edm::
InputTag>(
"muonCollectionTag")),
206 m_minTrackPt(
cfg.getParameter<double>(
"minTrackPt")),
207 m_maxTrackPt(
cfg.getParameter<double>(
"maxTrackPt")),
208 m_minTrackP(
cfg.getParameter<double>(
"minTrackP")),
209 m_maxTrackP(
cfg.getParameter<double>(
"maxTrackP")),
210 m_maxDxy(
cfg.getParameter<double>(
"maxDxy")),
211 m_minTrackerHits(
cfg.getParameter<
int>(
"minTrackerHits")),
212 m_maxTrackerRedChi2(
cfg.getParameter<double>(
"maxTrackerRedChi2")),
213 m_allowTIDTEC(
cfg.getParameter<
bool>(
"allowTIDTEC")),
214 m_minNCrossedChambers(
cfg.getParameter<
int>(
"minNCrossedChambers")),
215 m_minDT13Hits(
cfg.getParameter<
int>(
"minDT13Hits")),
216 m_minDT2Hits(
cfg.getParameter<
int>(
"minDT2Hits")),
217 m_minCSCHits(
cfg.getParameter<
int>(
"minCSCHits")),
218 m_writeTemporaryFile(
cfg.getParameter<
std::
string>(
"writeTemporaryFile")),
220 m_doAlignment(
cfg.getParameter<
bool>(
"doAlignment")),
221 m_strategy(
cfg.getParameter<
int>(
"strategy")),
222 m_residualsModel(
cfg.getParameter<
std::
string>(
"residualsModel")),
223 m_minAlignmentHits(
cfg.getParameter<
int>(
"minAlignmentHits")),
224 m_twoBin(
cfg.getParameter<
bool>(
"twoBin")),
225 m_combineME11(
cfg.getParameter<
bool>(
"combineME11")),
226 m_weightAlignment(
cfg.getParameter<
bool>(
"weightAlignment")),
227 m_reportFileName(
cfg.getParameter<
std::
string>(
"reportFileName")),
228 m_maxResSlopeY(
cfg.getParameter<double>(
"maxResSlopeY")),
229 m_createNtuple(
cfg.getParameter<
bool>(
"createNtuple")),
230 m_peakNSigma(
cfg.getParameter<double>(
"peakNSigma")),
231 m_BFieldCorrection(
cfg.getParameter<
int>(
"bFieldCorrection")),
232 m_doDT(
cfg.getParameter<
bool>(
"doDT")),
233 m_doCSC(
cfg.getParameter<
bool>(
"doCSC")),
234 m_useResiduals(
cfg.getParameter<
std::
string>(
"useResiduals")) {
279 m_ttree =
fs->make<TTree>(
"mual_ttree",
"mual_ttree");
305 return atoi(
s.c_str());
313 if (alignableMuon ==
nullptr)
314 throw cms::Exception(
"MuonAlignmentFromReference") <<
"doMuon must be set to True" << std::endl;
333 <<
"unrecognized residualsModel: \"" <<
m_residualsModel <<
"\"" << std::endl;
350 <<
"unrecognized useResiduals: \"" <<
m_useResiduals <<
"\"" << std::endl;
361 bool made_fitter =
false;
414 <<
"only DTChambers and CSCChambers can be aligned with this module" << std::endl;
420 int index = ali->geomDetId().rawId();
430 const auto& all_DT_chambers = alignableMuon->
DTChambers();
431 const auto& all_CSC_chambers = alignableMuon->
CSCChambers();
443 std::cout <<
"****** EVENT START *******" << std::endl;
455 std::cout <<
"JUST BEFORE LOOP OVER trajTrackPairs" << std::endl;
459 for (ConstTrajTrackPairCollection::const_iterator trajtrack = trajtracks.begin(); trajtrack != trajtracks.end();
473 std::cout <<
"JUST BEFORE muonResidualsFromTrack" << std::endl;
484 std::cout <<
"JUST AFTER muonResidualsFromTrack" << std::endl;
487 std::cout <<
"JUST BEFORE PROCESS" << std::endl;
490 std::cout <<
"JUST AFTER PROCESS" << std::endl;
495 std::cout <<
"JUST AFTER LOOP OVER trajTrackPairs" << std::endl;
540 std::vector<DetId> chamberIds = mrft.
chamberIds();
547 for (std::vector<DetId>::const_iterator chamberId = chamberIds.begin(); chamberId != chamberIds.end();
559 if (dt13 !=
nullptr && dt2 !=
nullptr) {
565 std::map<Alignable*, MuonResidualsTwoBin*>::const_iterator fitter =
592 std::cout <<
"processMuonResidualsFromTrack 6DOF dt13->residual() " << dt13->
residual()
610 fitter->second->fill(
charge, residdata);
625 if (dt13 !=
nullptr) {
630 std::map<Alignable*, MuonResidualsTwoBin*>::const_iterator fitter =
653 std::cout <<
"processMuonResidualsFromTrack 5DOF dt13->residual() " << dt13->
residual()
667 fitter->second->fill(
charge, residdata);
678 if (
csc !=
nullptr) {
688 std::map<Alignable*, MuonResidualsTwoBin*>::const_iterator fitter =
m_fitters.find(ali);
704 std::cout <<
"processMuonResidualsFromTrack 6DOFrphi csc->residual() " <<
csc->residual()
718 fitter->second->fill(
charge, residdata);
763 <<
"That's all!" << std::endl;
765 TStopwatch stop_watch;
772 std::cout <<
"readTmpFiles took " << stop_watch.CpuTime() <<
" sec" << std::endl;
781 std::cout <<
"selectResidualsPeaks took " << stop_watch.CpuTime() <<
" sec" << std::endl;
789 std::cout <<
"correctBField took " << stop_watch.CpuTime() <<
" sec" << std::endl;
798 std::cout <<
"fiducialCuts took " << stop_watch.CpuTime() <<
" sec" << std::endl;
807 std::cout <<
"fillNtuple took " << stop_watch.CpuTime() <<
" sec" << std::endl;
815 std::cout <<
"eraseNotSelectedResiduals took " << stop_watch.CpuTime() <<
" sec" << std::endl;
824 std::cout <<
"fitAndAlign took " << stop_watch.CpuTime() <<
" sec" << std::endl;
832 std::cout <<
"end: MuonAlignmentFromReference::terminate()" << std::endl;
845 report <<
"nan = None; NAN = None" << std::endl;
846 report <<
"nan = 0" << std::endl;
847 report <<
"reports = []" << std::endl;
848 report <<
"class ValErr:" << std::endl
849 <<
" def __init__(self, value, error, antisym):" << std::endl
850 <<
" self.value, self.error, self.antisym = value, error, antisym" << std::endl
852 <<
" def __repr__(self):" << std::endl
853 <<
" if self.antisym == 0.:" << std::endl
854 <<
" return \"%g +- %g\" % (self.value, self.error)" << std::endl
855 <<
" else:" << std::endl
856 <<
" return \"%g +- %g ~ %g\" % (self.value, self.error, self.antisym)" << std::endl
858 <<
"class Report:" << std::endl
859 <<
" def __init__(self, chamberId, postal_address, name):" << std::endl
860 <<
" self.chamberId, self.postal_address, self.name = chamberId, postal_address, name" << std::endl
861 <<
" self.status = \"NOFIT\"" << std::endl
862 <<
" self.fittype = None" << std::endl
864 <<
" def add_parameters(self, deltax, deltay, deltaz, deltaphix, deltaphiy, deltaphiz, loglikelihood, " 865 "numsegments, sumofweights, redchi2):" 867 <<
" self.status = \"PASS\"" << std::endl
868 <<
" self.deltax, self.deltay, self.deltaz, self.deltaphix, self.deltaphiy, self.deltaphiz = deltax, " 869 "deltay, deltaz, deltaphix, deltaphiy, deltaphiz" 871 <<
" self.loglikelihood, self.numsegments, self.sumofweights, self.redchi2 = loglikelihood, " 872 "numsegments, sumofweights, redchi2" 875 <<
" def add_stats(self, median_x, median_y, median_dxdz, median_dydz, mean30_x, mean30_y, mean20_dxdz, " 876 "mean50_dydz, mean15_x, mean15_y, mean10_dxdz, mean25_dydz, wmean30_x, wmean30_y, wmean20_dxdz, " 877 "wmean50_dydz, wmean15_x, wmean15_y, wmean10_dxdz, wmean25_dydz, stdev30_x, stdev30_y, stdev20_dxdz, " 878 "stdev50_dydz, stdev15_x, stdev15_y, stdev10_dxdz, stdev25_dydz):" 880 <<
" self.median_x, self.median_y, self.median_dxdz, self.median_dydz, self.mean30_x, self.mean30_y, " 881 "self.mean20_dxdz, self.mean50_dydz, self.mean15_x, self.mean15_y, self.mean10_dxdz, self.mean25_dydz, " 882 "self.wmean30_x, self.wmean30_y, self.wmean20_dxdz, self.wmean50_dydz, self.wmean15_x, self.wmean15_y, " 883 "self.wmean10_dxdz, self.wmean25_dydz, self.stdev30_x, self.stdev30_y, self.stdev20_dxdz, " 884 "self.stdev50_dydz, self.stdev15_x, self.stdev15_y, self.stdev10_dxdz, self.stdev25_dydz = median_x, " 885 "median_y, median_dxdz, median_dydz, mean30_x, mean30_y, mean20_dxdz, mean50_dydz, mean15_x, mean15_y, " 886 "mean10_dxdz, mean25_dydz, wmean30_x, wmean30_y, wmean20_dxdz, wmean50_dydz, wmean15_x, wmean15_y, " 887 "wmean10_dxdz, wmean25_dydz, stdev30_x, stdev30_y, stdev20_dxdz, stdev50_dydz, stdev15_x, stdev15_y, " 888 "stdev10_dxdz, stdev25_dydz" 891 <<
" def __repr__(self):" << std::endl
892 <<
" return \"<Report %s %s %s>\" % (self.postal_address[0], \" \".join(map(str, " 893 "self.postal_address[1:])), self.status)" 899 std::cout <<
"***** just after report.open" << std::endl;
903 std::cout <<
"***** Start loop over alignables" << std::endl;
905 std::vector<bool> selector = ali->alignmentParameters()->selector();
906 bool align_x = selector[0];
907 bool align_y = selector[1];
908 bool align_z = selector[2];
909 bool align_phix = selector[3];
910 bool align_phiy = selector[4];
911 bool align_phiz = selector[5];
912 int numParams = ((align_x ? 1 : 0) + (align_y ? 1 : 0) + (align_z ? 1 : 0) + (align_phix ? 1 : 0) +
913 (align_phiy ? 1 : 0) + (align_phiz ? 1 : 0));
916 std::vector<int> paramIndex;
917 int paramIndex_counter = -1;
919 paramIndex_counter++;
920 paramIndex.push_back(paramIndex_counter);
922 paramIndex_counter++;
923 paramIndex.push_back(paramIndex_counter);
925 paramIndex_counter++;
926 paramIndex.push_back(paramIndex_counter);
928 paramIndex_counter++;
929 paramIndex.push_back(paramIndex_counter);
931 paramIndex_counter++;
932 paramIndex.push_back(paramIndex_counter);
934 paramIndex_counter++;
935 paramIndex.push_back(paramIndex_counter);
937 DetId id = ali->geomDetId();
947 std::cout <<
"***** loop over alignables 1" << std::endl;
950 char wheel_label[][2] = {
"A",
"B",
"C",
"D",
"E"};
957 sprintf(cname,
"MBwh%sst%dsec%02d", wheel_label[chamberId.
wheel() + 2], chamberId.
station(), chamberId.
sector());
959 report <<
"reports.append(Report(" <<
id.rawId() <<
", (\"DT\", " << chamberId.
wheel() <<
", " 960 << chamberId.
station() <<
", " << chamberId.
sector() <<
"), \"" << cname <<
"\"))" << std::endl;
966 (chamberId.
endcap() == 1 ?
"p" :
"m"),
974 report <<
"reports.append(Report(" <<
id.rawId() <<
", (\"CSC\", " << chamberId.
endcap() <<
", " 975 << chamberId.
station() <<
", " << chamberId.
ring() <<
", " << chamberId.
chamber() <<
"), \"" << cname
976 <<
"\"))" << std::endl;
981 std::cout <<
"***** loop over alignables 2" << std::endl;
985 std::map<Alignable*, MuonResidualsTwoBin*>::const_iterator fitter =
m_fitters.find(thisali);
988 std::cout <<
"***** loop over alignables 3" << std::endl;
993 TStopwatch stop_watch;
998 std::cout <<
"=============================================================================================" 1001 std::cout <<
"Fitting " << cname << std::endl;
1004 report <<
"reports[-1].posNum = " << fitter->second->numResidualsPos() << std::endl;
1005 report <<
"reports[-1].negNum = " << fitter->second->numResidualsNeg() << std::endl;
1049 std::cout <<
"***** loop over alignables 4" << std::endl;
1056 std::cout <<
"***** loop over alignables 5" << std::endl;
1058 bool successful_fit = fitter->second->fit(thisali);
1061 std::cout <<
"***** loop over alignables 6 " << fitter->second->type() << std::endl;
1063 double loglikelihood = fitter->second->loglikelihood();
1064 double numsegments = fitter->second->numsegments();
1065 double sumofweights = fitter->second->sumofweights();
1066 double redchi2 = fitter->second->plot(cname, &rootDirectory, thisali);
1070 std::cout <<
"***** loop over alignables k5DOF" << std::endl;
1100 double gammaresid_value, gammaresid_error, gammaresid_antisym, gammaresslope_value, gammaresslope_error,
1101 gammaresslope_antisym;
1102 gammaresid_value = gammaresid_error = gammaresid_antisym = gammaresslope_value = gammaresslope_error =
1103 gammaresslope_antisym = 0.;
1118 report <<
"reports[-1].fittype = \"5DOF\"" << std::endl;
1119 report <<
"reports[-1].add_parameters(ValErr(" << deltax_value <<
", " << deltax_error <<
", " 1120 << deltax_antisym <<
"), \\" << std::endl
1121 <<
" None, \\" << std::endl
1122 <<
" ValErr(" << deltaz_value <<
", " << deltaz_error <<
", " 1123 << deltaz_antisym <<
"), \\" << std::endl
1124 <<
" ValErr(" << deltaphix_value <<
", " << deltaphix_error <<
", " 1125 << deltaphix_antisym <<
"), \\" << std::endl
1126 <<
" ValErr(" << deltaphiy_value <<
", " << deltaphiy_error <<
", " 1127 << deltaphiy_antisym <<
"), \\" << std::endl
1128 <<
" ValErr(" << deltaphiz_value <<
", " << deltaphiz_error <<
", " 1129 << deltaphiz_antisym <<
"), \\" << std::endl
1130 <<
" " << loglikelihood <<
", " << numsegments <<
", " << sumofweights
1131 <<
", " << redchi2 <<
")" << std::endl;
1132 report <<
"reports[-1].sigmaresid = ValErr(" << sigmaresid_value <<
", " << sigmaresid_error <<
", " 1133 << sigmaresid_antisym <<
")" << std::endl;
1134 report <<
"reports[-1].sigmaresslope = ValErr(" << sigmaresslope_value <<
", " << sigmaresslope_error
1135 <<
", " << sigmaresslope_antisym <<
")" << std::endl;
1139 report <<
"reports[-1].gammaresid = ValErr(" << gammaresid_value <<
", " << gammaresid_error <<
", " 1140 << gammaresid_antisym <<
")" << std::endl;
1141 report <<
"reports[-1].gammaresslope = ValErr(" << gammaresslope_value <<
", " << gammaresslope_error
1142 <<
", " << gammaresslope_antisym <<
")" << std::endl;
1167 <<
"None)" << std::endl;
1169 std::stringstream namesimple_x, namesimple_dxdz, nameweighted_x, nameweighted_dxdz;
1170 namesimple_x << cname <<
"_simple_x";
1171 namesimple_dxdz << cname <<
"_simple_dxdz";
1172 nameweighted_x << cname <<
"_weighted_x";
1173 nameweighted_dxdz << cname <<
"_weighted_dxdz";
1176 fitter->second->plotsimple(
1179 fitter->second->plotweighted(nameweighted_x.str(),
1184 fitter->second->plotweighted(nameweighted_dxdz.str(),
1191 if (successful_fit) {
1193 params[paramIndex[0]] = deltax_value;
1195 params[paramIndex[2]] = deltaz_value;
1197 params[paramIndex[3]] = deltaphix_value;
1199 params[paramIndex[4]] = deltaphiy_value;
1201 params[paramIndex[5]] = deltaphiz_value;
1207 std::cout <<
"***** loop over alignables k6DOF" << std::endl;
1249 double gammax_value, gammax_error, gammax_antisym, gammay_value, gammay_error, gammay_antisym,
1250 gammadxdz_value, gammadxdz_error, gammadxdz_antisym, gammadydz_value, gammadydz_error, gammadydz_antisym;
1251 gammax_value = gammax_error = gammax_antisym = gammay_value = gammay_error = gammay_antisym =
1252 gammadxdz_value = gammadxdz_error = gammadxdz_antisym = gammadydz_value = gammadydz_error =
1253 gammadydz_antisym = 0.;
1275 report <<
"reports[-1].fittype = \"6DOF\"" << std::endl;
1276 report <<
"reports[-1].add_parameters(ValErr(" << deltax_value <<
", " << deltax_error <<
", " 1277 << deltax_antisym <<
"), \\" << std::endl
1278 <<
" ValErr(" << deltay_value <<
", " << deltay_error <<
", " 1279 << deltay_antisym <<
"), \\" << std::endl
1280 <<
" ValErr(" << deltaz_value <<
", " << deltaz_error <<
", " 1281 << deltaz_antisym <<
"), \\" << std::endl
1282 <<
" ValErr(" << deltaphix_value <<
", " << deltaphix_error <<
", " 1283 << deltaphix_antisym <<
"), \\" << std::endl
1284 <<
" ValErr(" << deltaphiy_value <<
", " << deltaphiy_error <<
", " 1285 << deltaphiy_antisym <<
"), \\" << std::endl
1286 <<
" ValErr(" << deltaphiz_value <<
", " << deltaphiz_error <<
", " 1287 << deltaphiz_antisym <<
"), \\" << std::endl
1288 <<
" " << loglikelihood <<
", " << numsegments <<
", " << sumofweights
1289 <<
", " << redchi2 <<
")" << std::endl;
1290 report <<
"reports[-1].sigmax = ValErr(" << sigmax_value <<
", " << sigmax_error <<
", " << sigmax_antisym
1291 <<
")" << std::endl;
1292 report <<
"reports[-1].sigmay = ValErr(" << sigmay_value <<
", " << sigmay_error <<
", " << sigmay_antisym
1293 <<
")" << std::endl;
1294 report <<
"reports[-1].sigmadxdz = ValErr(" << sigmadxdz_value <<
", " << sigmadxdz_error <<
", " 1295 << sigmadxdz_antisym <<
")" << std::endl;
1296 report <<
"reports[-1].sigmadydz = ValErr(" << sigmadydz_value <<
", " << sigmadydz_error <<
", " 1297 << sigmadydz_antisym <<
")" << std::endl;
1301 report <<
"reports[-1].gammax = ValErr(" << gammax_value <<
", " << gammax_error <<
", " << gammax_antisym
1302 <<
")" << std::endl;
1303 report <<
"reports[-1].gammay = ValErr(" << gammay_value <<
", " << gammay_error <<
", " << gammay_antisym
1304 <<
")" << std::endl;
1305 report <<
"reports[-1].gammadxdz = ValErr(" << gammadxdz_value <<
", " << gammadxdz_error <<
", " 1306 << gammadxdz_antisym <<
")" << std::endl;
1307 report <<
"reports[-1].gammadydz = ValErr(" << gammadydz_value <<
", " << gammadydz_error <<
", " 1308 << gammadydz_antisym <<
")" << std::endl;
1347 std::stringstream namesimple_x, namesimple_y, namesimple_dxdz, namesimple_dydz, nameweighted_x,
1348 nameweighted_y, nameweighted_dxdz, nameweighted_dydz;
1349 namesimple_x << cname <<
"_simple_x";
1350 namesimple_y << cname <<
"_simple_y";
1351 namesimple_dxdz << cname <<
"_simple_dxdz";
1352 namesimple_dydz << cname <<
"_simple_dydz";
1353 nameweighted_x << cname <<
"_weighted_x";
1354 nameweighted_y << cname <<
"_weighted_y";
1355 nameweighted_dxdz << cname <<
"_weighted_dxdz";
1356 nameweighted_dydz << cname <<
"_weighted_dydz";
1360 fitter->second->plotsimple(
1362 fitter->second->plotsimple(
1365 fitter->second->plotweighted(nameweighted_x.str(),
1370 fitter->second->plotweighted(nameweighted_y.str(),
1375 fitter->second->plotweighted(nameweighted_dxdz.str(),
1380 fitter->second->plotweighted(nameweighted_dydz.str(),
1387 if (successful_fit) {
1389 params[paramIndex[0]] = deltax_value;
1391 params[paramIndex[1]] = deltay_value;
1393 params[paramIndex[2]] = deltaz_value;
1395 params[paramIndex[3]] = deltaphix_value;
1397 params[paramIndex[4]] = deltaphiy_value;
1399 params[paramIndex[5]] = deltaphiz_value;
1405 std::cout <<
"***** loop over alignables k6DOFrphi" << std::endl;
1439 double gammaresid_value, gammaresid_error, gammaresid_antisym, gammaresslope_value, gammaresslope_error,
1440 gammaresslope_antisym;
1441 gammaresid_value = gammaresid_error = gammaresid_antisym = gammaresslope_value = gammaresslope_error =
1442 gammaresslope_antisym = 0.;
1456 report <<
"reports[-1].fittype = \"6DOFrphi\"" << std::endl;
1457 report <<
"reports[-1].add_parameters(ValErr(" << deltax_value <<
", " << deltax_error <<
", " 1458 << deltax_antisym <<
"), \\" << std::endl
1459 <<
" ValErr(" << deltay_value <<
", " << deltay_error <<
", " 1460 << deltay_antisym <<
"), \\" << std::endl
1461 <<
" ValErr(" << deltaz_value <<
", " << deltaz_error <<
", " 1462 << deltaz_antisym <<
"), \\" << std::endl
1463 <<
" ValErr(" << deltaphix_value <<
", " << deltaphix_error <<
", " 1464 << deltaphix_antisym <<
"), \\" << std::endl
1465 <<
" ValErr(" << deltaphiy_value <<
", " << deltaphiy_error <<
", " 1466 << deltaphiy_antisym <<
"), \\" << std::endl
1467 <<
" ValErr(" << deltaphiz_value <<
", " << deltaphiz_error <<
", " 1468 << deltaphiz_antisym <<
"), \\" << std::endl
1469 <<
" " << loglikelihood <<
", " << numsegments <<
", " << sumofweights
1470 <<
", " << redchi2 <<
")" << std::endl;
1471 report <<
"reports[-1].sigmaresid = ValErr(" << sigmaresid_value <<
", " << sigmaresid_error <<
", " 1472 << sigmaresid_antisym <<
")" << std::endl;
1473 report <<
"reports[-1].sigmaresslope = ValErr(" << sigmaresslope_value <<
", " << sigmaresslope_error
1474 <<
", " << sigmaresslope_antisym <<
")" << std::endl;
1478 report <<
"reports[-1].gammaresid = ValErr(" << gammaresid_value <<
", " << gammaresid_error <<
", " 1479 << gammaresid_antisym <<
")" << std::endl;
1480 report <<
"reports[-1].gammaresslope = ValErr(" << gammaresslope_value <<
", " << gammaresslope_error
1481 <<
", " << gammaresslope_antisym <<
")" << std::endl;
1491 << fitter->second->wmean(
1495 << fitter->second->wmean(
1499 << fitter->second->wmean(
1503 << fitter->second->wmean(
1510 <<
"None)" << std::endl;
1512 std::stringstream namesimple_x, namesimple_dxdz, nameweighted_x, nameweighted_dxdz;
1513 namesimple_x << cname <<
"_simple_x";
1514 namesimple_dxdz << cname <<
"_simple_dxdz";
1515 nameweighted_x << cname <<
"_weighted_x";
1516 nameweighted_dxdz << cname <<
"_weighted_dxdz";
1519 fitter->second->plotsimple(
1522 fitter->second->plotweighted(nameweighted_x.str(),
1527 fitter->second->plotweighted(nameweighted_dxdz.str(),
1534 if (successful_fit) {
1536 params[paramIndex[0]] = deltax_value;
1538 params[paramIndex[1]] = deltay_value;
1540 params[paramIndex[2]] = deltaz_value;
1542 params[paramIndex[3]] = deltaphix_value;
1544 params[paramIndex[4]] = deltaphiy_value;
1546 params[paramIndex[5]] = deltaphiz_value;
1550 if (successful_fit) {
1552 oneortwo.push_back(ali);
1554 oneortwo.push_back(thisali);
1558 std::cout <<
"MINUIT fit failed!" << std::endl;
1560 report <<
"reports[-1].status = \"MINUITFAIL\"" << std::endl;
1563 for (
int i = 0;
i < numParams;
i++)
1567 oneortwo.push_back(ali);
1569 oneortwo.push_back(thisali);
1574 std::cout <<
"Too few hits!" << std::endl;
1576 report <<
"reports[-1].status = \"TOOFEWHITS\"" << std::endl;
1579 for (
int i = 0;
i < numParams;
i++)
1583 oneortwo.push_back(ali);
1585 oneortwo.push_back(thisali);
1590 ali->setAlignmentParameters(parnew);
1592 ali->alignmentParameters()->setValid(
true);
1595 std::cout << cname <<
" fittime= " << stop_watch.CpuTime() <<
" sec" << std::endl;
1613 file = fopen((*fileName).c_str(),
"r");
1614 if (
file ==
nullptr)
1616 <<
"file \"" << *
fileName <<
"\" can't be opened (doesn't exist?)" << std::endl;
1618 fread(&
size,
sizeof(
int), 1,
file);
1622 <<
" fitters (probably corresponds to the wrong alignment job)" << std::endl;
1627 unsigned int index_toread;
1628 fread(&index_toread,
sizeof(
unsigned int), 1,
file);
1629 if (*
index != index_toread)
1631 <<
"file \"" << *
fileName <<
"\" has index " << index_toread <<
" at position " <<
i 1632 <<
", but this job is expecting " << *
index <<
" (probably corresponds to the wrong alignment job)" 1645 fwrite(&
size,
sizeof(
int), 1,
file);
1650 unsigned int index_towrite = *
index;
1651 fwrite(&index_towrite,
sizeof(
unsigned int), 1,
file);
1695 int vars_index[10] = {0, 1};
1761 sprintf(cname,
"MB%+d/%d/%02d", chamberId.
wheel(), chamberId.
station(), chamberId.
sector());
1766 (chamberId.
endcap() == 1 ?
"+" :
"-"),
1802 for (; residual != fitter->
residualsPos_end(); ++residual, ++residual_ok) {
1840 std::map<Alignable*, bool> already_seen;
1843 bool parsing_error =
false;
1848 parsing_error =
true;
1850 if (!parsing_error &&
barrel) {
1864 parsing_error =
true;
1867 bool wheel_digit =
false;
1877 parsing_error =
true;
1880 parsing_error =
true;
1884 bool station_digit =
false;
1888 station_digit =
true;
1892 parsing_error =
true;
1895 parsing_error =
true;
1899 bool sector_digit =
false;
1903 sector_digit =
true;
1907 parsing_error =
true;
1909 if (!parsing_error) {
1910 bool no_such_chamber =
false;
1912 if (wheel < -2 || wheel > 2)
1913 no_such_chamber =
true;
1914 if (station < 1 || station > 4)
1915 no_such_chamber =
true;
1916 if (
station == 4 && (sector < 1 || sector > 14))
1917 no_such_chamber =
true;
1918 if (
station < 4 && (sector < 1 || sector > 12))
1919 no_such_chamber =
true;
1921 if (no_such_chamber)
1923 <<
"reference chamber doesn't exist: " << (*name) << std::endl;
1926 for (
const auto& ali : all_DT_chambers) {
1927 if (ali->geomDetId().rawId() ==
id.rawId()) {
1928 std::map<Alignable*, bool>::const_iterator
trial = already_seen.find(ali);
1929 if (
trial == already_seen.end()) {
1931 already_seen[ali] =
true;
1938 if (!parsing_error &&
endcap) {
1952 parsing_error =
true;
1955 bool station_digit =
false;
1959 station_digit =
true;
1965 parsing_error =
true;
1968 parsing_error =
true;
1972 bool ring_digit =
false;
1980 parsing_error =
true;
1983 parsing_error =
true;
1987 bool chamber_digit =
false;
1991 chamber_digit =
true;
1995 parsing_error =
true;
1997 if (!parsing_error) {
1998 bool no_such_chamber =
false;
2002 if (station < 1 || station > 4)
2003 no_such_chamber =
true;
2004 if (
station == 1 && (ring < 1 || ring > 4))
2005 no_such_chamber =
true;
2006 if (
station > 1 && (ring < 1 || ring > 2))
2007 no_such_chamber =
true;
2008 if (
station == 1 && (chamber < 1 || chamber > 36))
2009 no_such_chamber =
true;
2010 if (
station > 1 &&
ring == 1 && (chamber < 1 || chamber > 18))
2011 no_such_chamber =
true;
2012 if (
station > 1 &&
ring == 2 && (chamber < 1 || chamber > 36))
2013 no_such_chamber =
true;
2015 if (no_such_chamber)
2017 <<
"reference chamber doesn't exist: " << (*name) << std::endl;
2020 for (
const auto& ali : all_CSC_chambers) {
2021 if (ali->geomDetId().rawId() ==
id.rawId()) {
2022 std::map<Alignable*, bool>::const_iterator
trial = already_seen.find(ali);
2023 if (
trial == already_seen.end()) {
2025 already_seen[ali] =
true;
2034 <<
"reference chamber name is malformed: " << (*name) << std::endl;
void initialize(const edm::EventSetup &iSetup, AlignableTracker *alignableTracker, AlignableMuon *alignableMuon, AlignableExtras *extras, AlignmentParameterStore *alignmentParameterStore) override
Call at beginning of job (must be implemented in derived class)
edm::InputTag m_muonCollectionTag
int station() const
Return the station number.
double m_maxTrackerRedChi2
long m_counter_trackerchi2
int trackerNumHits() const
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > m_MagFieldToken
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
long m_counter_trackertidtec
const edm::ESGetToken< DetIdAssociator, DetIdAssociatorRecord > m_DetIdToken
align::Alignables CSCChambers()
align::Alignables m_alignables
void setAlignmentPositionError(const align::Alignables &alivec, double valshift, double valrot)
Set Alignment position error.
long m_counter_station123aligning
int m_minNCrossedChambers
void selectResidualsPeaks()
long m_counter_station123dt13hits
const align::Alignables & alignables(void) const
get all alignables
void applyParameters(void)
Obsolete: Use AlignableNavigator::alignableDetFromDetId and alignableFromAlignableDet.
std::string m_reportFileName
~MuonAlignmentFromReference() override
std::vector< bool >::const_iterator residualsPos_ok_begin() const
Interface/Base class for alignment algorithms, each alignment algorithm has to be derived from this c...
constexpr Detector det() const
get the detector field from this detid
long m_counter_trackmomentum
define event information passed to algorithms
std::string m_residualsModel
AlignmentParameterStore * m_alignmentParameterStore
void processMuonResidualsFromTrack(MuonResidualsFromTrack &mrft)
int number(std::string s)
MuonResidualsFitter::MuonAlignmentTreeRow m_tree_row
double pt() const
track transverse momentum
std::vector< unsigned int > m_indexes
int charge() const
track electric charge
std::string chamberPrettyNameFromId(unsigned int idx)
align::Alignables DTChambers()
const reco::Track * getTrack()
const MuonResidualsFromTrack::BuilderToken m_builderToken
void startNewLoop() override
AlignableNavigator * m_alignableNavigator
long m_counter_cscaligning
void eraseNotSelectedResiduals()
Abs< T >::type abs(const T &t)
double normalizedChi2() const
const edm::ESGetToken< Propagator, TrackingComponentsRecord > m_propToken
void parseReference(align::Alignables &reference, const align::Alignables &all_DT_chambers, const align::Alignables &all_CSC_chambers)
void write(FILE *file, int which=0)
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector's numbering enum) ...
MuonAlignmentFromReference(const edm::ParameterSet &cfg, edm::ConsumesCollector &iC)
const edm::ESGetToken< CSCGeometry, MuonGeometryRecord > m_cscGeometryToken
std::map< unsigned int, MuonResidualsTwoBin * > m_fitterOrder
long m_counter_station123valid
std::vector< std::string > m_reference
std::map< Alignable *, MuonResidualsTwoBin * > m_fitters
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
long m_counter_station4hits
void run(const edm::EventSetup &iSetup, const EventInfo &eventInfo) override
Run the algorithm (must be implemented in derived class)
CLHEP::HepVector AlgebraicVector
long m_counter_station123
TFileDirectory mkdir(const std::string &dir, const std::string &descr="")
create a new subdirectory
std::vector< double * >::const_iterator residualsPos_begin() const
void eraseNotSelectedResiduals()
const std::vector< DetId > chamberIds() const
long m_counter_station123dt2hits
const DetId & geomDetId() const
constexpr uint32_t rawId() const
get the raw id
long m_counter_totchambers
void read(FILE *file, int which=0)
double pz() const
z coordinate of momentum vector
std::string m_writeTemporaryFile
std::string m_useResiduals
long m_counter_station4aligning
long m_counter_station4valid
std::vector< Alignable * > Alignables
std::map< Alignable *, Alignable * > m_me11map
const edm::ESGetToken< GlobalTrackingGeometry, GlobalTrackingGeometryRecord > m_globTackingToken
std::vector< double * >::const_iterator residualsPos_end() const
int wheel() const
Return the wheel number.
virtual void terminate()
Called at end of job (must be implemented in derived class)
CLHEP::HepSymMatrix AlgebraicSymMatrix
std::vector< std::string > m_readTemporaryFiles
bool numeric(std::string s)
AlignableDetOrUnitPtr chamberAlignable() const
#define DEFINE_EDM_PLUGIN(factory, type, name)
bool contains_TIDTEC() const
eventInfo
add run, event number and lumi section
long m_counter_trackerhits
void selectPeakResiduals(double nsigma, int nvar, int *vars)
virtual AlignmentParameters * cloneFromSelected(const AlgebraicVector &par, const AlgebraicSymMatrix &cov) const =0
Constructor of the full muon geometry.
long m_counter_minchambers
MuonChamberResidual * chamberResidual(DetId chamberId, int type)
std::vector< ConstTrajTrackPair > ConstTrajTrackPairCollection