MY RESEARCH LOG

8 Mar 2012

Wow !! Happy New Year.... logging into this page after real long time. While doing work with python it helps if we know some stuff.....
  • variable : a variable is a name that refers to a location in computer memory where some value is stored. variables are usually created in an assignment statement like name = "stupido" . variable names can be arbitrarily long, must begin with a letter, can contain numbers, letters or even an underscore sign. a programmer can choose any name for his variable but for some reserved names called keywords. variables have an associated "type" which determines the possible values it can hold. usually the first assignment decides the type of a variable, which is fixed for the rest of program. for example
           name = 'stupid' # create a variable by assigment to a string value.  "name" can henceforth hold string values only
           name = 'silly' # legal statement
           name = 5 # illegal: and Integer value assigned to a variable of type string
          
    Python has five standard variable types:
    • Numbers
    • String
    • Tuple
    • Dictionary

  • statement : statement is a unit of code executable by the python interpreter. statement can either be a print statement which displays something on screen or it can be an assignment statement which either creates and give values to a new variable or simply assign a new value to an existing variable. an assignment statement leads to no output on screen.
  • script : a sequence of the statements constitute a script.
  • operators : Untill I find a better definition, operators for me just mean the special symbols in python that represent the computations e.g +, -. * and /.
  • operands : the values or variables on which the the operators act are called operands.
  • expression : an expression is formed by combining all the chips described above i.e values, variables and operators. extending this to extreme even a single value or a single variable are legal expression.

5 Sep 2011

I have been very busy whole last week with the generator level comparison of MC generated using SHERPA and MADGRAPH, and will try to summarize the attempt in a later post. Right now, its about getting on with the installation, usage and mastering of python-qt4 package. Here is how we install it on ubuntu:
sudo apt-get install python-qt4
After installation, now is the time to test it. I copy pasted following code into a file called "firstProgram.py"
import sys
from PyQt4.QtGui import *
application = QApplication(sys.argv)
button = QPushButton("Around the world !!!", None)
button.show()
application.exec_()
Running it creates a small window... containing the text. Here's a demonstrative project for a newbie like my self.
  • Project1: Create a widget with a title, a picture icon, a tooltip text and a push-button for clean exit.
#!/usr/bin/python
#project.py
import sys
from PyQt4 import QtGui
from PyQt4 import QtCore
class project(QtGui.QWidget):
    def __init__(self,parent=None):
        QtGui.QWidget.__init__(self,parent)
        self.setGeometry(300,300,250,150)
        self.setWindowTitle('MyProjectTitle')
        self.setWindowIcon(QtGui.QIcon('/home/anil/globe.png'))
        self.setToolTip('My first project page')
        quit = QtGui.QPushButton('Close',self)
        quit.setGeometry(10,10,60,35)
        self.connect(quit, QtCore.SIGNAL('clicked()'),
                     QtGui.qApp,QtCore.SLOT('quit()'))
app = QtGui.QApplication(sys.argv)
proj = project()
proj.show()
sys.exit(app.exec_())

  • Project2: Code a widget which prompts a message upon exiting, and ask you if you were sure about quitting.
import sys
from PyQt4 import QtGui
from PyQt4 import QtCore

class winWdMsgBox(QtGui.QWidget):
    def __init__(self,parent=None):
        QtGui.QWidget.__init__(self, parent)
        self.setGeometry(300,300,250,150)
        self.setWindowTitle("ExitWithMessage")
                
    def closeEvent(self, event):
        reply = QtGui.QMessageBox.question(
            self, 'Message',
            "Are you sure to quit?", QtGui.QMessageBox.Yes | 
            QtGui.QMessageBox.No, QtGui.QMessageBox.No)
       if reply == QtGui.QMessageBox.Yes:
            event.accept()
        else:
            event.ignore()
apl = QtGui.QApplication(sys.argv)
mb = winWdMsgBox()
mb.show()
sys.exit(apl.exec_())

Both these are motivated by very good tutorials on the internet, and both these work fine.

27 August 2011

Some root stuff: The fast fourier transform or FFT package can be installed and configured as follows on root
sudo apt-get install gfortran ncurses-dev libpcre3-dev xlibmesa-glu-dev libglew1.5-dev libftgl-dev libmysqlclient-dev libfftw3-dev cfitsio-dev graphviz-dev libavahi-compat-libdnssd-dev libldap-dev python-dev libxml2-dev libssl-dev libgsl0-dev gfortran-multilib gfortran-doc gfortran-4.4-multilib gfortran-4.4-doc libgfortran3-dbg glew-utils ocaml-findlib libtool-doc gsl-bin build-essential libjpeg62-dev libtiff4-dev  libgnutls-dev libgmp3-dev libxmu-dev libpng12-dev libkrb5-dev
freeglut3-devlibxmu-dev libgif-dev libgif-dev libiodbc2 libiodbc2-dev
./configure linux --enable-minuit2 --enable-roofit --enable-table --enable-gdml --enable-pgsql --enable-mysql --enable-builtin-zlib --enable-builtin-pcre --enable-builtin-freetype --enable-builtin-glew
The pileup studies seem like too tough an issue to begin with. I today tried to understand the trigInfoSummary producer module in CMSSW. It is a small code written by Mike Hidreth and produces the objects of type <std::vector< PileupSummaryInfo> >. The produce method begins by picking up the true pile-up information from the mixing module using a tag "mix" whose object type is PileupMixingContent. If the information is present and valid then,
  • Obtain the vector of bunch crossing ids.
  • Obtain a parallel vector of number of interactions.
If "mixing truth" information is not available, we work with the SimVertices. Now the module pick objects of type CrossingFrame using the tag "mix" + simHitLabel. Then the simVertex collection is used to create mix collection of type: std::auto_ptr<MixCollection >( new MixCollection==(cfSimVertexes.product()) ). Now we loop over this collection and store the bunch-crossing id, and store the number of interactions per crossing. Afterwards the code tries to collect the particles that look like as being in the tracking volume, it picks up the objects of type "TrackingParticleCollection" and "TrackingVertexCollection".

26 August 2011

The news today: Lovedeep tells that it is possible to have more than one type of pat-jet-collections (i.e separate patJets for iterativeCone, antiKT or KT algos) in pattuples. They key is to use "addJetCollection()" functionality. It appears cool and can be used when we are going to compare different type of Jet algorithms. DelPanj tree should also be optimized to be capable of storing more than one possible pat jet collections in branches. Requires some more changes to the code architecture...i am afraid.

25 August 2011

The big job for now is validation of the good old DelPanj tree. Trying to run PAT on CMSSW_4_2_0 RelVal samples which can then be fed into the tree maker code. The idea to begin with is
cmsrel CMSSW_4_2_5
cd CMSSW_4_2_5/src/
cvs co -r CMSSW_4_2_5 PhysicsTools/PatAlgos/test/patTuple_standard_cfg.py
cd PhysicsTools/PatAlgos/test/
cmsenv
cmsRun patTuple_standard_cfg.py
###Store the file in castor for future use.
rfcp pat_Output_File.root /castor/cern.ch/user/a/anil/
I have run the tree maker over the pat samples. All branches look fine except the ones carrying the un-corrected information about the jets. This is how I tried it:
double pt = -99, px = -99, py = -99, pz = -99, en = -99;
//Check if uncorrected jet information is available.
bool isAvailable = jet->jecSetAvailable("Uncorrected");
if(isAvailable){
LorentzVector uncorrP4 = jet->correctedP4("Uncorrected");
pt = uncorrP4.pt();
px = uncorrP4.px();
py = uncorrP4.py();
pz = uncorrP4.pz();
en = uncorrP4.e();
}
JetUnCorrPt_.push_back(pt);
JetUnCorrPx_.push_back(px);
JetUnCorrPy_.push_back(py);
JetUnCorrPz_.push_back(pz);
JetUnCorrEnt_.push_back(en);
But each of these branches are filled with "-99" i.e. default value.

22 August 2011

It's been almost a week since we landed at Madison. Most of the setting are done. Presented recently at one of the Mjj meetings, the talk was prepared in hurry and I would like to revisit the plots for better understanding. So now two tasks at hand, *PileUp Studies *V+Jets Effort for 2011 and Mjj I have spent a few hours poring over the literature and found some good twiki pages: The third is an effort by the Higgs people to look at the lepton selection and mass reconstruction in presence of high pileup.

2 August 2011

Here is the code which solved my missing et problems:
  using namespace edm;
  edm::Handle<HepMCProduct> event;
  iEvent.getByLabel("generator", event);
  const HepMC::GenEvent *myGenEvent = event->GetEvent();
  std::vector<const HepMC::GenParticle*> stableGenColl; 
  std::vector<const HepMC::GenParticle*> stableLepColl; 
  std::vector<const HepMC::GenParticle*> allNeutrinoColl;
  TLorentzVector momsum(0., 0., 0., 0.);  
  TLorentzVector neutrino_momsum(0., 0., 0., 0.);
 
  for(HepMC::GenEvent::particle_const_iterator giter = myGenEvent->particles_begin(); giter != myGenEvent->particles_end(); ++giter) {
    if((*giter)->status()!=1) continue;
    nsigned int pid = fabs((*giter)->pdg_id());   
    stableGenColl.push_back(*giter);
    bool isneutrino = false;
    TLorentzVector mom((*giter)->momentum().x(), (*giter)->momentum().y(), (*giter)->momentum().z(), (*giter)->momentum().t());
    if(pid==12 || pid ==14||pid==16){
      isneutrino = true;
      allNeutrinoColl.push_back(*giter);
      neutrino_momsum += mom;
     }
    if(!(isneutrino)){
      momsum += mom;
      double pt = (*giter)->momentum().perp();
      double eta = (*giter)->momentum().eta();
      if(pid==11 || pid ==13||pid==15){
   if(!(eta>3 ||pt<20))
     stableLepColl.push_back(*giter);
      }
    }
  }
  //Fill Missing Et stuff.
  if(!allNeutrinoColl.size()){std::cout<<"retiring"<<std::endl ;return;}
  hCounter_->Fill(1);
  TLorentzVector met(-1*momsum.Px(),-1*momsum.Py(),0,momsum.Pt());
  TLorentzVector neutrino(allNeutrinoColl[0]->momentum().x(), allNeutrinoColl[0]->momentum().y(), 
           allNeutrinoColl[0]->momentum().z(), allNeutrinoColl[0]->momentum().t()); 

  TLorentzVector neutrinomet(neutrino_momsum.Px(),neutrino_momsum.Py(),0,neutrino_momsum.Pt());
  TLorentzVector recmet(met);
  hRecMet_->Fill(recmet.Et());
  hNeutrinoPt_->Fill(neutrino.Et());
  hNeutrinoSumPt_->Fill(neutrino_momsum.Pt());
  hNeutrinoMet_->Fill(neutrinomet.Et());

1 August 2011

Working with official Madgraph W+Jet(s) samples. A strange problem, dont know if it is a problem at all. I loop over all the particles in the HepMCProduct, and calculate the resultant momentum over all the stable particles, except the neutrinos. We call transverse component of this grand resultant : GenMet and assume that for W+Jet MC, the GenMet must be same as the neutrino momentum. And since this is the generator level, I dont see any reason why they should not be exactly same with in the computer accuracy.... but my first plots are unfortunately not converging with this hypothesis. The GenMet distribution looks more diffused, and has a greater width compared to neutrino momentum. I posted this problem along with the plots (click) and got following reply
Hi Anil,
    There could be two sources that can produce this effect:
- the presence of other neutrinos in the event, from decays
- the fact that you sum particles within a large but not infinite eta range (most likely -5,+5)
I believe if you try to enlarge the eta range the agreement will become better and better.
Let me know
Cheers
Giulio
This is how i plan to check the first point, we can fill the number of neutrinos in the event. The difference between the genMET and the JetPt must be strongly correlated to the numNeutrinos in the event. Also I am recalculating the missing Et, but this time only those neutrinoes would be kept out of calculation of grand transverse momentum who have W boson as mother. Doing so should make the agreement better if Giulio's first argument is true. About the eta range, I have removed all explicit eta cuts while calculating the grand momentum... so assuming the HepMCProduct does not have any implicit cut on eta...., the second argument does not apply to our case.

Well it turns out that first point is correct. If take care of the fact that not all the neutrinos come from W decay, the agreement is perfect. Here (click) is the plot. Only those neutrinoes have been excluded from the sumPT calculation who come from the W bosons. Although the genMET now agrees with the NeutrinoPt, but the basic philosopy of calculation appears wrong to me at this moment. Instead of calculating genMET while using the non-W neutrinoes as if they are like other particles, we should not use any neutrino information. Rather we must calculate the genMET as I was doing originally but compare it to the resultant transverse momentum of all the neutrinoes in the event and not just the one coming from W decay. Here is the picture for that (click)

Finally here is the summary: GenMet must be compared to resultant of transverse momentum of all the neutrinoes in the event and not just the one which comes from the W decay. This is shown in this plot (click).

31 July 2011

I had a severe problem with the Njet distribution. It appears that we have forgotten to switch on the GenJet reconstruction and so the multiplicities seen here are merely the handiwork of the coding.... Since each of our generation has the MCProduct stored, I think that simply running GenJet reconstruction over them will do the job. Currently I am moving to the Madgraph samples generated centrally by CMS and placed in the Wisconsin Tier2 center, while lovedeep resolves this issue.

14 July 2011

Beginning of the day and I am already bogged down by TFileService (most traitorous):
@@@@ Checking shared library for missing symbols: libanalysisMjjAnalysis.so
tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so: 
                                undefined reference to `TFileDirectory::_cd(std::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, bool) const'
tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so: 
                                undefined reference to `TH1AddDirectorySentry::~TH1AddDirectorySentry()'
tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so: 
                                undefined reference to `TH1AddDirectorySentry::TH1AddDirectorySentry()'
tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so: 
                                undefined reference to `typeinfo for TFileService'
collect2: ld returned 1 exit status
gmake: *** [tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so] Error 1
Still have to figure out what worked for me, but I got rid of the problem by updating my BuildFile.xml blindly to following state:
<use name="FWCore/Framework"/>
<use name="FWCore/PluginManager"/>
<use name="FWCore/ParameterSet"/>
<use name="hepmc"/>
<use name="SimDataFormats/GeneratorProducts"/>
<use name="clhep"/>
<use name="root"/>
<use name="DataFormats/Candidate">
<use name="FWCore/Framework"/>
<use name="FWCore/PluginManager"/>
<use name="FWCore/ParameterSet"/>
<use name="boost"/>
<use name="rootcore"/>
<use name="rootrflx"/>
<use name="PhysicsTools/SelectorUtils"/>
<use name="PhysicsTools/FWLite"/>
<use name="DataFormats/Math"/>
<use name="DataFormats/Candidate"/>
<use name="CommonTools/UtilAlgos"/>
<flags EDM_PLUGIN="1"/>
<export>
   <lib name="1"/>
</export>

Still I get following exception at runtime:

%MSG-s CMSException:  AfterModConstruction 14-Jul-2011 11:27:52 CEST pre-events
cms::Exception caught in cmsRun
---- Configuration BEGIN
Error occurred while creating for module of type 'MjjAnalysis' with label 'demo'
---- NotFound BEGIN
Service Request unable to find requested service with compiler type name ' 12TFileService'.
---- NotFound END
---- Configuration END
Well it is about the way you call TFileService in your configuration. They correct way is:
process.TFileService = cms.Service("TFileService",
      fileName = cms.string("histo.root"),
      closeFileFast = cms.untracked.bool(True)
  )

I am in habit of forgetting the bussines of sorting the vector of TLorentzVectors: Here is the recipie->

//Define the predicate class
class PtGreater {
public:
  template <typename T> bool operator () (const T& i, const T& j) {  return (i.Pt() > j.Pt());  }
};

void foo(){
std::vector<TLorentzVector> someVector;
//Fill the vector
---------------
----------------
//sort the vector
std::sort(someVector.begin(),someVector.end(),PtGreater());
}

13 July 2011

Today I started the day trying to fetch the generation information using the gen particle collection in the followign manner:
   edm::Handle<reco::GenParticleCollection> genParticleHandle;
   iEvent.getByLabel(genPartLabel_, genParticleHandle)
   const reco::GenParticleCollection* genColl= &(*genParticleHandle);
   reco::GenParticleCollection::const_iterator geni = genColl->begin();
   for(; geni!=genColl->end();geni++){
   //---------------------   
   }
so far it is fine, but as soon as we try to put following line into the loop,
   for(; geni!=genColl->end();geni++){
    reco::GenParticle gen = *geni;
   }
the scram build fails with following error report:
**** ERROR: Failed parsing file: src/analysis/MjjAnalysis/BuildFile.xml.
>> Local Products Rules ..... started
>> Local Products Rules ..... done
>> Entering Package analysis/MjjAnalysis
>>Creating project symlinks
Entering library rule at analysis/MjjAnalysis
>> Compiling /afs/cern.ch/user/a/anil/scratch0/Gen/CMSSW_4_1_7/src/analysis/MjjAnalysis/src/MjjAnalysis.cc 
>> Building shared library tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so
@@@@ Checking shared library for missing symbols: libanalysisMjjAnalysis.so
tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so: undefined reference to `reco::GenParticle::~GenParticle()'
tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so: undefined reference to `vtable for reco::GenParticle'
collect2: ld returned 1 exit status
gmake: *** [tmp/slc5_amd64_gcc434/src/analysis/MjjAnalysis/src/analysisMjjAnalysis/libanalysisMjjAnalysis.so] Error 1
And to be further sure here is the buildfile that I am using:
<use name="FWCore/Framework"/>
<use name="FWCore/PluginManager"/>
<use name="FWCore/ParameterSet"/>
<use name="hepmc"/>
<use name="SimDataFormats/GeneratorProducts"/>
<use name="clhep"/>
<use name="root"/>
<use name="DataFormats/Candidate">
<flags EDM_PLUGIN="1"/>
<export>
   <lib name="1"/>
</export>
Although I could not find the reason for the crash but have found a way out instead.... use HepMCProduct instead of gen particle collection. Here is the algorithm:
  • Gather the HepMCProduct information. edm::Handle event; iEvent.getByLabel(hepMcColl_, evt);
  • Get the gen event, const HepMC::GenEvent *myGenEvent = evt->GetEvent();

  • We have a particle iterator which can be used as follows:
          HepMC::GenEvent::particle_const_iterator iter;
          iter = myGenEvent->particles_begin();
          //======== some stuff with iter
          iter++;
     

Here are the line that worked for me:

 using namespace edm;
   edm::Handle<HepMCProduct> event;
   iEvent.getByLabel(hepMcColl_, event);
   const HepMC::GenEvent *myGenEvent = event->GetEvent();
   std::vector<const HepMC::GenParticle*> stableGenColl;

   for(HepMC::GenEvent::particle_const_iterator giter = myGenEvent->particles_begin(); giter != myGenEvent->particles_end(); ++giter) {
    if((*giter)->status()!=1)continue;
    stableGenColl.push_back(*giter);
   }
and ofcourse there was an include:
#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"

11 July 2011

Well the first problem is the crowding of the space due to the LSFJOB directories. Well here is the plan to deal with it. I write a script that runs in background and looks for the LSFJOB directories and them kicks the ones with STDOUT and other files into the temp area. But before that I submit a new set of jobs.

Edit | Attach | Watch | Print version | History: r28 < r27 < r26 < r25 < r24 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r28 - 2012-03-09 - AnilPratapSingh
 
    • Cern Search Icon Cern Search
    • TWiki Search Icon TWiki Search
    • Google Search Icon Google Search

    Main All webs login

This site is powered by the TWiki collaboration platform Powered by PerlCopyright & 2008-2019 by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
Ideas, requests, problems regarding TWiki? Send feedback