CMS 3D CMS Logo

/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/Alignment/MillePedeAlignmentAlgorithm/src/Mille.cc

Go to the documentation of this file.
00001 
00010 #include "Mille.h"
00011 
00012 #include <iostream>
00013 
00014 //___________________________________________________________________________
00015 
00016 Mille::Mille(const char *outFileName, bool asBinary, bool writeZero) : 
00017   myOutFile(outFileName, (asBinary ? (std::ios::binary | std::ios::out) : std::ios::out)),
00018   myAsBinary(asBinary), myWriteZero(writeZero), myBufferPos(-1), myHasSpecial(false)
00019 {
00020   // opens outFileName, by default as binary file
00021 
00022   // Instead myBufferPos(-1), myHasSpecial(false) and the following two lines
00023   // we could call newSet() and kill()...
00024   myBufferInt[0]   = 0;
00025   myBufferFloat[0] = 0.;
00026 
00027   if (!myOutFile.is_open()) {
00028     std::cerr << "Mille::Mille: Could not open " << outFileName 
00029               << " as output file." << std::endl;
00030   }
00031 }
00032 
00033 //___________________________________________________________________________
00034 
00035 Mille::~Mille()
00036 {
00037   // closes file
00038   myOutFile.close();
00039 }
00040 
00041 //___________________________________________________________________________
00042 
00043 void Mille::mille(int NLC, const float *derLc,
00044                   int NGL, const float *derGl, const int *label,
00045                   float rMeas, float sigma)
00046 {
00047   if (sigma <= 0.) return;
00048   if (myBufferPos == -1) this->newSet(); // start, e.g. new track
00049   if (!this->checkBufferSize(NLC, NGL)) return;
00050 
00051   // first store measurement
00052   ++myBufferPos;
00053   myBufferFloat[myBufferPos] = rMeas;
00054   myBufferInt  [myBufferPos] = 0;
00055 
00056   // store local derivatives and local 'lables' 1,...,NLC
00057   for (int i = 0; i < NLC; ++i) {
00058     if (derLc[i] || myWriteZero) { // by default store only non-zero derivatives
00059       ++myBufferPos;
00060       myBufferFloat[myBufferPos] = derLc[i]; // local derivatives
00061       myBufferInt  [myBufferPos] = i+1;      // index of local parameter
00062     }
00063   }
00064 
00065   // store uncertainty of measurement in between locals and globals
00066   ++myBufferPos;
00067   myBufferFloat[myBufferPos] = sigma;
00068   myBufferInt  [myBufferPos] = 0;
00069 
00070   // store global derivatives and their lables
00071   for (int i = 0; i < NGL; ++i) {
00072     if (derGl[i] || myWriteZero) { // by default store only non-zero derivatives
00073       if ((label[i] > 0 || myWriteZero) && label[i] <= myMaxLabel) { // and for valid labels
00074         ++myBufferPos;
00075         myBufferFloat[myBufferPos] = derGl[i]; // global derivatives
00076         myBufferInt  [myBufferPos] = label[i]; // index of global parameter
00077       } else {
00078         std::cerr << "Mille::mille: Invalid label " << label[i] 
00079                   << " <= 0 or > " << myMaxLabel << std::endl; 
00080       }
00081     }
00082   }
00083 }
00084 
00085 //___________________________________________________________________________
00086 void Mille::special(int nSpecial, const float *floatings, const int *integers)
00087 {
00088   if (nSpecial == 0) return;
00089   if (myBufferPos == -1) this->newSet(); // start, e.g. new track
00090   if (myHasSpecial) {
00091     std::cerr << "Mille::special: Special values already stored for this record."
00092               << std::endl; 
00093     return;
00094   }
00095   if (!this->checkBufferSize(nSpecial, 0)) return;
00096   myHasSpecial = true; // after newSet() (Note: MILLSP sets to buffer position...)
00097 
00098   //  myBufferFloat[.]  | myBufferInt[.]
00099   // ------------------------------------
00100   //      0.0           |      0
00101   //  -float(nSpecial)  |      0
00102   //  The above indicates special data, following are nSpecial floating and nSpecial integer data.
00103 
00104   ++myBufferPos; // zero pair
00105   myBufferFloat[myBufferPos] = 0.;
00106   myBufferInt  [myBufferPos] = 0;
00107 
00108   ++myBufferPos; // nSpecial and zero
00109   myBufferFloat[myBufferPos] = -nSpecial; // automatic conversion to float
00110   myBufferInt  [myBufferPos] = 0;
00111 
00112   for (int i = 0; i < nSpecial; ++i) {
00113     ++myBufferPos;
00114     myBufferFloat[myBufferPos] = floatings[i];
00115     myBufferInt  [myBufferPos] = integers[i];
00116   }
00117 }
00118 
00119 //___________________________________________________________________________
00120 
00121 void Mille::kill()
00122 {
00123   // reset buffers, i.e. kill derivatives accumulated for current set
00124   myBufferPos = -1;
00125 }
00126 
00127 //___________________________________________________________________________
00128 
00129 void Mille::end()
00130 {
00131   // write set of derivatives with same local parameters to file
00132   if (myBufferPos > 0) { // only if anything stored...
00133     const int numWordsToWrite = (myBufferPos + 1)*2;
00134 
00135     if (myAsBinary) {
00136       myOutFile.write(reinterpret_cast<const char*>(&numWordsToWrite), 
00137                       sizeof(numWordsToWrite));
00138       myOutFile.write(reinterpret_cast<char*>(myBufferFloat), 
00139                       (myBufferPos+1) * sizeof(myBufferFloat[0]));
00140       myOutFile.write(reinterpret_cast<char*>(myBufferInt), 
00141                       (myBufferPos+1) * sizeof(myBufferInt[0]));
00142     } else {
00143       myOutFile << numWordsToWrite << "\n";
00144       for (int i = 0; i < myBufferPos+1; ++i) {
00145         myOutFile << myBufferFloat[i] << " ";
00146       }
00147       myOutFile << "\n";
00148       
00149       for (int i = 0; i < myBufferPos+1; ++i) {
00150         myOutFile << myBufferInt[i] << " ";
00151       }
00152       myOutFile << "\n";
00153     }
00154   }
00155   myBufferPos = -1; // reset buffer for next set of derivatives
00156 }
00157 
00158 //___________________________________________________________________________
00159 
00160 void Mille::newSet()
00161 {
00162   // initilise for new set of locals, e.g. new track
00163   myBufferPos = 0;
00164   myHasSpecial = false;
00165   myBufferFloat[0] = 0.0;
00166   myBufferInt  [0] = 0;   // position 0 used as error counter
00167 }
00168 
00169 //___________________________________________________________________________
00170 
00171 bool Mille::checkBufferSize(int nLocal, int nGlobal)
00172 {
00173   // enough space for next nLocal + nGlobal derivatives incl. measurement?
00174 
00175   if (myBufferPos + nLocal + nGlobal + 2 >= myBufferSize) {
00176     ++(myBufferInt[0]); // increase error count
00177     std::cerr << "Mille::checkBufferSize: Buffer too short (" 
00178               << myBufferSize << "),"
00179               << "\n need space for nLocal (" << nLocal<< ")"
00180               << "/nGlobal (" << nGlobal << ") local/global derivatives, " 
00181               << myBufferPos + 1 << " already stored!"
00182               << std::endl;
00183     return false;
00184   } else {
00185     return true;
00186   }
00187 }