CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch12/src/CondCore/EcalPlugins/plugins/EcalPedestalsPyWrapper.cc

Go to the documentation of this file.
00001 #include "CondFormats/EcalObjects/interface/EcalPedestals.h"
00002 #include "CondTools/Ecal/interface/EcalPedestalsXMLTranslator.h"
00003 #include "CondTools/Ecal/interface/EcalCondHeader.h"
00004 #include "TH2F.h"
00005 #include "TCanvas.h"
00006 #include "TStyle.h"
00007 #include "TLine.h"
00008 #include "DataFormats/EcalDetId/interface/EBDetId.h"
00009 #include "DataFormats/EcalDetId/interface/EEDetId.h"
00010 
00011 #include "CondCore/Utilities/interface/PayLoadInspector.h"
00012 #include "CondCore/Utilities/interface/InspectorPythonWrapper.h"
00013 
00014 #include <string>
00015 #include <sstream>
00016 #include <algorithm>
00017 #include <numeric>
00018 #include <iterator>
00019 #include <boost/ref.hpp>
00020 #include <boost/bind.hpp>
00021 #include <boost/function.hpp>
00022 #include <boost/iterator/transform_iterator.hpp>
00023 
00024 #include <fstream>
00025 
00026 namespace {
00027   struct Printer {
00028     void doit(EcalPedestal const & item) {
00029       for (int i=1; i<4; i++)
00030         ss << item.mean(i) << "," << item.rms(i) <<";";
00031       ss << " ";
00032     }
00033     std::stringstream ss;
00034   };
00035 }
00036 
00037 namespace cond {
00038 
00039   namespace ecalped {
00040     enum Quantity { mean_x12=1, mean_x6=2, mean_x3=3 };
00041     enum How { singleChannel, bySuperModule, all};
00042 
00043     float average(EcalPedestals const & peds, Quantity q) {
00044       return std::accumulate(
00045                              boost::make_transform_iterator(peds.barrelItems().begin(),bind(&EcalPedestal::mean,_1,q)),
00046                              boost::make_transform_iterator(peds.barrelItems().end(),bind(&EcalPedestal::mean,_1,q)),
00047                              0.)/float(peds.barrelItems().size());
00048     }
00049 
00050     void extractAverage(EcalPedestals const & peds, Quantity q, std::vector<int> const &,  std::vector<float> & result) {
00051       result.resize(1);
00052       result[0] = average(peds,q);
00053     }
00054     
00055     void extractSuperModules(EcalPedestals const & peds, Quantity q, std::vector<int> const & which,  std::vector<float> & result) {
00056       // bho...
00057     }
00058 
00059     void extractSingleChannel(EcalPedestals const & peds, Quantity q, std::vector<int> const & which,  std::vector<float> & result) {
00060       for (unsigned int i=0; i<which.size();i++) {
00061         // absolutely arbitraty
00062         if ((unsigned int) (which[i])<  peds.barrelItems().size())
00063           result.push_back( peds.barrelItems()[which[i]].mean(q));
00064       }
00065     }
00066 
00067         typedef boost::function<void(EcalPedestals const & peds, Quantity q, std::vector<int> const & which,  std::vector<float> & result)> PedExtractor;
00068   }
00069 
00070   template<>
00071   struct ExtractWhat<EcalPedestals> {
00072 
00073     ecalped::Quantity m_quantity;
00074     ecalped::How m_how;
00075     std::vector<int> m_which;
00076 
00077     ecalped::Quantity const & quantity() const { return m_quantity;}
00078     ecalped::How const & how() const { return m_how;}
00079     std::vector<int> const & which() const { return m_which;}
00080  
00081 
00082     void set_quantity( ecalped::Quantity i) { m_quantity=i;}
00083     void set_how(ecalped::How i) {m_how=i;}
00084     void set_which(std::vector<int> & i) { m_which.swap(i);}
00085   };
00086 
00087 
00088   template<>
00089   class ValueExtractor<EcalPedestals>: public  BaseValueExtractor<EcalPedestals> {
00090   public:
00091 
00092     static ecalped::PedExtractor & extractor(ecalped::How how) {
00093       static  ecalped::PedExtractor fun[3] = { 
00094         ecalped::PedExtractor(ecalped::extractSingleChannel),
00095         ecalped::PedExtractor(ecalped::extractSuperModules),
00096         ecalped::PedExtractor(ecalped::extractAverage)
00097               };
00098       return fun[how];
00099     }
00100 
00101     typedef EcalPedestals Class;
00102     typedef ExtractWhat<Class> What;
00103     static What what() { return What();}
00104 
00105     ValueExtractor(){}
00106     ValueExtractor(What const & what)
00107       : m_what(what)
00108     {
00109       // here one can make stuff really complicated... (select mean rms, 12,6,1)
00110       // ask to make average on selected channels...
00111     }
00112 
00113     void compute(Class const & it){
00114       std::vector<float> res;
00115       extractor(m_what.how())(it,m_what.quantity(),m_what.which(),res);
00116       swap(res);
00117     }
00118 
00119   private:
00120     What  m_what;
00121     
00122   };
00123 
00124 
00125   template<>
00126   std::string
00127   PayLoadInspector<EcalPedestals>::dump() const {
00128     std::stringstream ss;
00129     EcalCondHeader header;
00130     ss<<EcalPedestalsXMLTranslator::dumpXML(header,object());
00131     return ss.str();
00132   }  // dump
00133 
00134   template<>
00135   std::string PayLoadInspector<EcalPedestals>::summary() const {
00136     std::cout << "***************************************"<< std::endl;
00137     std::stringstream ss;
00138     ss << "sizes="
00139        << object().barrelItems().size() <<","
00140        << object().endcapItems().size() <<";";
00141     ss << std::endl;
00142     return ss.str();
00143   }  // summary
00144 
00145 
00146   // return the real name of the file including extension...
00147   template<>
00148   std::string PayLoadInspector<EcalPedestals>::plot(std::string const & filename,
00149                                                     std::string const &, 
00150                                                     std::vector<int> const&, 
00151                                                     std::vector<float> const& ) const {
00152     gStyle->SetPalette(1);
00153     //    TCanvas canvas("CC map","CC map",840,600);
00154     TCanvas canvas("CC map","CC map",800,1200);
00155  
00156     float xmi[3] = {0.0 , 0.22, 0.78};
00157     float xma[3] = {0.22, 0.78, 1.00};
00158     TPad*** pad = new TPad**[6];
00159     for (int gId = 0; gId < 6; gId++) {
00160       pad[gId] = new TPad*[3];
00161       for (int obj = 0; obj < 3; obj++) {
00162         float yma = 1.- (0.17 * gId);
00163         float ymi = yma - 0.15;
00164         pad[gId][obj] = new TPad(Form("p_%i_%i", obj, gId),Form("p_%i_%i", obj, gId),
00165                                  xmi[obj], ymi, xma[obj], yma);
00166         pad[gId][obj]->Draw();
00167       }
00168     }
00169 
00170     const int kGains       = 3;
00171     const int gainValues[3] = {12, 6, 1};
00172     const int kSides       = 2;
00173     const int kBarlRings   = EBDetId::MAX_IETA;
00174     const int kBarlWedges  = EBDetId::MAX_IPHI;
00175     const int kEndcWedgesX = EEDetId::IX_MAX;
00176     const int kEndcWedgesY = EEDetId::IY_MAX;
00177 
00178     TH2F** barrel_m = new TH2F*[3];
00179     TH2F** endc_p_m = new TH2F*[3];
00180     TH2F** endc_m_m = new TH2F*[3];
00181     TH2F** barrel_r = new TH2F*[3];
00182     TH2F** endc_p_r = new TH2F*[3];
00183     TH2F** endc_m_r = new TH2F*[3];
00184     for (int gainId = 0; gainId < kGains; gainId++) {
00185       barrel_m[gainId] = new TH2F(Form("EBm%i",gainId),Form("mean %i EB",gainValues[gainId]),360,0,360, 170, -85,85);
00186       endc_p_m[gainId] = new TH2F(Form("EE+m%i",gainId),Form("mean %i EE+",gainValues[gainId]),100,1,101,100,1,101);
00187       endc_m_m[gainId] = new TH2F(Form("EE-m%i",gainId),Form("mean %i EE-",gainValues[gainId]),100,1,101,100,1,101);
00188       barrel_r[gainId] = new TH2F(Form("EBr%i",gainId),Form("rms %i EB",gainValues[gainId]),360,0,360, 170, -85,85);
00189       endc_p_r[gainId] = new TH2F(Form("EE+r%i",gainId),Form("rms %i EE+",gainValues[gainId]),100,1,101,100,1,101);
00190       endc_m_r[gainId] = new TH2F(Form("EE-r%i",gainId),Form("rms %i EE-",gainValues[gainId]),100,1,101,100,1,101);
00191     }
00192 
00193     for (int sign=0; sign < kSides; sign++) {
00194       int thesign = sign==1 ? 1:-1;
00195 
00196       for (int ieta=0; ieta<kBarlRings; ieta++) {
00197         for (int iphi=0; iphi<kBarlWedges; iphi++) {
00198           EBDetId id((ieta+1)*thesign, iphi+1);
00199           float y = -1 - ieta;
00200           if(sign == 1) y = ieta;
00201           barrel_m[0]->Fill(iphi, y, object()[id.rawId()].mean_x12);
00202           barrel_r[0]->Fill(iphi, y, object()[id.rawId()].rms_x12);
00203           barrel_m[1]->Fill(iphi, y, object()[id.rawId()].mean_x6);
00204           barrel_r[1]->Fill(iphi, y, object()[id.rawId()].rms_x6);
00205           barrel_m[2]->Fill(iphi, y, object()[id.rawId()].mean_x1);
00206           barrel_r[2]->Fill(iphi, y, object()[id.rawId()].rms_x1);
00207         }  // iphi
00208       }   // ieta
00209 
00210       for (int ix=0; ix<kEndcWedgesX; ix++) {
00211         for (int iy=0; iy<kEndcWedgesY; iy++) {
00212           if (! EEDetId::validDetId(ix+1,iy+1,thesign)) continue;
00213           EEDetId id(ix+1,iy+1,thesign);
00214           if (thesign==1) {
00215             endc_p_m[0]->Fill(ix+1,iy+1,object()[id.rawId()].mean_x12);
00216             endc_p_r[0]->Fill(ix+1,iy+1,object()[id.rawId()].rms_x12);
00217             endc_p_m[1]->Fill(ix+1,iy+1,object()[id.rawId()].mean_x6);
00218             endc_p_r[1]->Fill(ix+1,iy+1,object()[id.rawId()].rms_x6);
00219             endc_p_m[2]->Fill(ix+1,iy+1,object()[id.rawId()].mean_x1);
00220             endc_p_r[2]->Fill(ix+1,iy+1,object()[id.rawId()].rms_x1);
00221           }
00222           else{ 
00223             endc_m_m[0]->Fill(ix+1,iy+1,object()[id.rawId()].mean_x12);
00224             endc_m_r[0]->Fill(ix+1,iy+1,object()[id.rawId()].rms_x12);
00225             endc_m_m[1]->Fill(ix+1,iy+1,object()[id.rawId()].mean_x6);
00226             endc_m_r[1]->Fill(ix+1,iy+1,object()[id.rawId()].rms_x6);
00227             endc_m_m[2]->Fill(ix+1,iy+1,object()[id.rawId()].mean_x1);
00228             endc_m_r[2]->Fill(ix+1,iy+1,object()[id.rawId()].rms_x1);
00229           }
00230         }  // iy
00231       }   // ix
00232     }    // side
00233 
00234     //canvas.cd(1);
00235     float bmin[3] ={0.7, 0.5, 0.4};
00236     float bmax[3] ={1.7, 1.0, 0.8};
00237     float emin[3] ={1.5, 0.8, 0.4};
00238     float emax[3] ={2.5, 1.5, 0.8};
00239     TLine* l = new TLine(0., 0., 0., 0.);
00240     l->SetLineWidth(1);
00241     int ixSectorsEE[202] = {
00242       62, 62, 61, 61, 60, 60, 59, 59, 58, 58, 56, 56, 46, 46, 44, 44, 43, 43, 42, 42, 
00243       41, 41, 40, 40, 41, 41, 42, 42, 43, 43, 44, 44, 46, 46, 56, 56, 58, 58, 59, 59, 
00244       60, 60, 61, 61, 62, 62,  0,101,101, 98, 98, 96, 96, 93, 93, 88, 88, 86, 86, 81, 
00245       81, 76, 76, 66, 66, 61, 61, 41, 41, 36, 36, 26, 26, 21, 21, 16, 16, 14, 14,  9,
00246       9,  6,  6,  4,  4,  1,  1,  4,  4,  6,  6,  9,  9, 14, 14, 16, 16, 21, 21, 26, 
00247       26, 36, 36, 41, 41, 61, 61, 66, 66, 76, 76, 81, 81, 86, 86, 88, 88, 93, 93, 96, 
00248       96, 98, 98,101,101,  0, 62, 66, 66, 71, 71, 81, 81, 91, 91, 93,  0, 62, 66, 66, 
00249       91, 91, 98,  0, 58, 61, 61, 66, 66, 71, 71, 76, 76, 81, 81,  0, 51, 51,  0, 44, 
00250       41, 41, 36, 36, 31, 31, 26, 26, 21, 21,  0, 40, 36, 36, 11, 11,  4,  0, 40, 36, 
00251       36, 31, 31, 21, 21, 11, 11,  9,  0, 46, 46, 41, 41, 36, 36,  0, 56, 56, 61, 61, 66, 66};
00252 
00253     int iySectorsEE[202] = {
00254       51, 56, 56, 58, 58, 59, 59, 60, 60, 61, 61, 62, 62, 61, 61, 60, 60, 59, 59, 58, 
00255       58, 56, 56, 46, 46, 44, 44, 43, 43, 42, 42, 41, 41, 40, 40, 41, 41, 42, 42, 43, 
00256       43, 44, 44, 46, 46, 51,  0, 51, 61, 61, 66, 66, 76, 76, 81, 81, 86, 86, 88, 88, 
00257       93, 93, 96, 96, 98, 98,101,101, 98, 98, 96, 96, 93, 93, 88, 88, 86, 86, 81, 81, 
00258       76, 76, 66, 66, 61, 61, 41, 41, 36, 36, 26, 26, 21, 21, 16, 16, 14, 14,  9,  9, 
00259       6,  6,  4,  4,  1,  1,  4,  4,  6,  6,  9,  9, 14, 14, 16, 16, 21, 21, 26, 26, 
00260       36, 36, 41, 41, 51,  0, 46, 46, 41, 41, 36, 36, 31, 31, 26, 26,  0, 51, 51, 56, 
00261       56, 61, 61,  0, 61, 61, 66, 66, 71, 71, 76, 76, 86, 86, 88,  0, 62,101,  0, 61, 
00262       61, 66, 66, 71, 71, 76, 76, 86, 86, 88,  0, 51, 51, 56, 56, 61, 61,  0, 46, 46, 
00263       41, 41, 36, 36, 31, 31, 26, 26,  0, 40, 31, 31, 16, 16,  6,  0, 40, 31, 31, 16, 16,  6};
00264   
00265     for (int gId = 0; gId < 3; gId++) {
00266       pad[gId][0]->cd();
00267       endc_m_m[gId]->SetStats(0);
00268       endc_m_m[gId]->SetMaximum(225);
00269       endc_m_m[gId]->SetMinimum(175);
00270       endc_m_m[gId]->Draw("colz");
00271       for ( int i=0; i<201; i=i+1) {
00272         if ( (ixSectorsEE[i]!=0 || iySectorsEE[i]!=0) && 
00273              (ixSectorsEE[i+1]!=0 || iySectorsEE[i+1]!=0) ) {
00274           l->DrawLine(ixSectorsEE[i], iySectorsEE[i], 
00275                       ixSectorsEE[i+1], iySectorsEE[i+1]);
00276           l->SetLineWidth(0.2);
00277         }
00278       }
00279       pad[gId + 3][0]->cd();
00280       endc_m_r[gId]->SetStats(0);
00281       endc_m_r[gId]->SetMaximum(emax[gId]);
00282       endc_m_r[gId]->SetMinimum(emin[gId]);
00283       endc_m_r[gId]->Draw("colz");
00284       for ( int i=0; i<201; i=i+1) {
00285         if ( (ixSectorsEE[i]!=0 || iySectorsEE[i]!=0) && 
00286              (ixSectorsEE[i+1]!=0 || iySectorsEE[i+1]!=0) ) {
00287           l->DrawLine(ixSectorsEE[i], iySectorsEE[i], 
00288                       ixSectorsEE[i+1], iySectorsEE[i+1]);
00289         }
00290       }
00291     //canvas.cd(2);
00292       pad[gId][1]->cd();
00293       barrel_m[gId]->SetStats(0);
00294       barrel_m[gId]->SetMaximum(225);
00295       barrel_m[gId]->SetMinimum(175);
00296       barrel_m[gId]->Draw("colz");
00297       for(int i = 0; i <17; i++) {
00298         Double_t x = 20.+ (i *20);
00299         l = new TLine(x,-85.,x,86.);
00300         l->Draw();
00301       }
00302       l = new TLine(0.,0.,360.,0.);
00303       l->Draw();
00304       pad[gId + 3][1]->cd();
00305       barrel_r[gId]->SetStats(0);
00306       barrel_r[gId]->SetMaximum(bmax[gId]);
00307       barrel_r[gId]->SetMinimum(bmin[gId]);
00308       barrel_r[gId]->Draw("colz");
00309       for(int i = 0; i <17; i++) {
00310         Double_t x = 20.+ (i *20);
00311         l = new TLine(x,-85.,x,86.);
00312         l->Draw();
00313       }
00314       l = new TLine(0.,0.,360.,0.);
00315       l->Draw();
00316     //canvas.cd(3);
00317       pad[gId][2]->cd();
00318       endc_p_m[gId]->SetStats(0);
00319       endc_p_m[gId]->SetMaximum(225);
00320       endc_p_m[gId]->SetMinimum(175);
00321       endc_p_m[gId]->Draw("colz");
00322       for ( int i=0; i<201; i=i+1) {
00323         if ( (ixSectorsEE[i]!=0 || iySectorsEE[i]!=0) && 
00324              (ixSectorsEE[i+1]!=0 || iySectorsEE[i+1]!=0) ) {
00325           l->DrawLine(ixSectorsEE[i], iySectorsEE[i], 
00326                       ixSectorsEE[i+1], iySectorsEE[i+1]);
00327         }
00328       }
00329       pad[gId + 3][2]->cd();
00330       endc_p_r[gId]->SetStats(0);
00331       endc_p_r[gId]->SetMaximum(emax[gId]);
00332       endc_p_r[gId]->SetMinimum(emin[gId]);
00333       endc_p_r[gId]->Draw("colz");
00334       for ( int i=0; i<201; i=i+1) {
00335         if ( (ixSectorsEE[i]!=0 || iySectorsEE[i]!=0) && 
00336              (ixSectorsEE[i+1]!=0 || iySectorsEE[i+1]!=0) ) {
00337           l->DrawLine(ixSectorsEE[i], iySectorsEE[i], 
00338                       ixSectorsEE[i+1], iySectorsEE[i+1]);
00339         }
00340       }
00341     }
00342 
00343     canvas.SaveAs(filename.c_str());
00344     return filename;
00345   }  // plot
00346 }
00347 
00348 namespace condPython {
00349   template<>
00350   void defineWhat<EcalPedestals>() {
00351     using namespace boost::python;
00352     enum_<cond::ecalped::Quantity>("Quantity")
00353       .value("mean_x12",cond::ecalped::mean_x12)
00354       .value("mean_x6",  cond::ecalped::mean_x6)
00355       .value("mean_x3", cond::ecalped::mean_x3)
00356       ;
00357     enum_<cond::ecalped::How>("How")
00358       .value("singleChannel",cond::ecalped::singleChannel)
00359       .value("bySuperModule",cond::ecalped::bySuperModule) 
00360       .value("all",cond::ecalped::all)
00361       ;
00362 
00363     typedef cond::ExtractWhat<EcalPedestals> What;
00364     class_<What>("What",init<>())
00365       .def("set_quantity",&What::set_quantity)
00366       .def("set_how",&What::set_how)
00367       .def("set_which",&What::set_which)
00368       .def("quantity",&What::quantity, return_value_policy<copy_const_reference>())
00369       .def("how",&What::how, return_value_policy<copy_const_reference>())
00370       .def("which",&What::which, return_value_policy<copy_const_reference>())
00371       ;
00372   }
00373 }
00374 
00375 
00376 
00377 PYTHON_WRAPPER(EcalPedestals,EcalPedestals);