File indexing completed on 2023-10-25 09:56:55
0001
0002
0003
0004
0005
0006 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0007 #include "FWCore/Framework/interface/EventSetup.h"
0008 #include "FWCore/Framework/interface/ESHandle.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011
0012 #include "MagneticField/VolumeBasedEngine/interface/MagGeometry.h"
0013 #include "MagneticField/GeomBuilder/test/stubs/MagGeometryExerciser.h"
0014 #include "MagneticField/Engine/interface/MagneticField.h"
0015 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0016 #include "MagneticField/GeomBuilder/src/MagGeoBuilderFromDDD.h"
0017 #include "MagneticField/VolumeBasedEngine/interface/VolumeBasedMagneticField.h"
0018
0019 #include <iostream>
0020 #include <vector>
0021
0022 class testMagGeometryAnalyzer : public edm::one::EDAnalyzer<> {
0023 public:
0024
0025 testMagGeometryAnalyzer(const edm::ParameterSet& pset);
0026
0027
0028 ~testMagGeometryAnalyzer() override = default;
0029
0030
0031 void analyze(const edm::Event& event, const edm::EventSetup& eventSetup) override;
0032
0033 void endJob() override {}
0034
0035 private:
0036 void testGrids(const std::vector<MagVolume6Faces const*>& bvol, const VolumeBasedMagneticField* field);
0037
0038 const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magfieldToken_;
0039 };
0040
0041 testMagGeometryAnalyzer::testMagGeometryAnalyzer(const edm::ParameterSet&)
0042 : magfieldToken_(esConsumes<MagneticField, IdealMagneticFieldRecord>()) {}
0043
0044 void testMagGeometryAnalyzer::analyze(const edm::Event& event, const edm::EventSetup& eventSetup) {
0045 const edm::ESHandle<MagneticField>& magfield = eventSetup.getHandle(magfieldToken_);
0046
0047 const VolumeBasedMagneticField* field = dynamic_cast<const VolumeBasedMagneticField*>(magfield.product());
0048 const MagGeometry* geom = field->field;
0049
0050
0051
0052 MagGeometryExerciser exe(geom);
0053 exe.testFindVolume(1000000);
0054
0055
0056
0057
0058
0059
0060
0061 if (true) {
0062 edm::LogVerbatim("MagGeometry") << "***TEST GRIDS: barrel volumes: " << geom->barrelVolumes().size();
0063 testGrids(geom->barrelVolumes(), field);
0064
0065 edm::LogVerbatim("MagGeometry") << "***TEST GRIDS: endcap volumes: " << geom->endcapVolumes().size();
0066 testGrids(geom->endcapVolumes(), field);
0067 }
0068 }
0069
0070 #include "MagneticField/VolumeGeometry/interface/MagVolume6Faces.h"
0071 #include "VolumeGridTester.h"
0072
0073 void testMagGeometryAnalyzer::testGrids(const std::vector<MagVolume6Faces const*>& bvol,
0074 const VolumeBasedMagneticField* field) {
0075 static std::map<std::string, int> nameCalls;
0076
0077 for (std::vector<MagVolume6Faces const*>::const_iterator i = bvol.begin(); i != bvol.end(); i++) {
0078 if ((*i)->copyno != 1) {
0079 continue;
0080 }
0081
0082 const MagProviderInterpol* prov = (**i).provider();
0083 if (prov == 0) {
0084 edm::LogVerbatim("MagGeometry") << (*i)->volumeNo << " No interpolator; skipping ";
0085 continue;
0086 }
0087 VolumeGridTester tester(*i, prov, field);
0088 if (tester.testInside())
0089 edm::LogVerbatim("MagGeometry") << "testGrids: success: " << (**i).volumeNo;
0090 else
0091 edm::LogVerbatim("MagGeometry") << "testGrids: ERROR: " << (**i).volumeNo;
0092 }
0093 }
0094
0095 #include "FWCore/Framework/interface/MakerMacros.h"
0096 DEFINE_FWK_MODULE(testMagGeometryAnalyzer);