Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-19 01:44:22

0001 /*
0002  * Authors: William Nash (original), Sven Dildick (adapted)
0003  */
0004 
0005 #include "L1Trigger/CSCTriggerPrimitives/interface/CSCPatternBank.h"
0006 
0007 #include <TCanvas.h>
0008 #include <TGraph.h>
0009 #include <TGraphErrors.h>
0010 #include <TF1.h>
0011 #include <TSystem.h>
0012 #include <TH1F.h>
0013 #include <TFile.h>
0014 #include <TString.h>
0015 
0016 //c++
0017 #include <fstream>
0018 #include <iostream>
0019 #include <iomanip>
0020 #include <math.h>
0021 #include <memory>
0022 #include <bitset>
0023 
0024 using namespace std;
0025 
0026 // uncertainty per layer in half strips 1/sqrt(12)
0027 const float HALF_STRIP_ERROR = 0.288675;
0028 const unsigned halfpatternwidth = (CSCConstants::CLCT_PATTERN_WIDTH - 1) / 2;
0029 
0030 //labels of all envelopes
0031 const int DEBUG = 0;
0032 
0033 // forward declarations
0034 class CSCPattern {
0035 public:
0036   CSCPattern(unsigned int id, const CSCPatternBank::LCTPattern& pat);
0037   ~CSCPattern() {}
0038 
0039   void printCode(const unsigned code) const;
0040   void getLayerPattern(const unsigned code, unsigned layerPattern[CSCConstants::NUM_LAYERS]) const;
0041   int recoverPatternCCCombination(const unsigned code,
0042                                   int code_hits[CSCConstants::NUM_LAYERS][CSCConstants::CLCT_PATTERN_WIDTH]) const;
0043   string getName() const { return name_; }
0044   unsigned getId() const { return id_; }
0045 
0046 private:
0047   const unsigned int id_;
0048   string name_;
0049   CSCPatternBank::LCTPattern pat_;
0050 };
0051 
0052 CSCPattern::CSCPattern(unsigned int id, const CSCPatternBank::LCTPattern& pat)
0053     : id_(id), name_(std::to_string(id)), pat_(pat) {}
0054 
0055 //given a code, prints out how it looks within the pattern
0056 void CSCPattern::printCode(const unsigned code) const {
0057   // comparator code per layer
0058   unsigned layerPattern[CSCConstants::NUM_LAYERS];
0059   getLayerPattern(code, layerPattern);
0060 
0061   std::cout << "Pattern " << id_ << ", Code " << code << " " << std::bitset<12>(code) << std::endl;
0062   for (unsigned int j = 0; j < CSCConstants::NUM_LAYERS; j++) {
0063     unsigned trueCounter = 0;  //for each layer, should only have 3
0064     std::cout << "L" << j + 1 << ": ";
0065     for (unsigned int i = 0; i < CSCConstants::CLCT_PATTERN_WIDTH; i++) {
0066       if (!pat_[j][i]) {
0067         printf("-");
0068       } else {
0069         trueCounter++;
0070         if (trueCounter == layerPattern[j])
0071           printf("X");
0072         else
0073           printf("_");
0074       }
0075     }
0076     std::cout << " -> " << std::bitset<2>(layerPattern[j]) << std::endl;
0077   }
0078 
0079   printf("    ");
0080   for (unsigned int i = 0; i < CSCConstants::CLCT_PATTERN_WIDTH; i++) {
0081     printf("%1x", i);
0082   }
0083   printf("\n");
0084 }
0085 
0086 void CSCPattern::getLayerPattern(const unsigned code, unsigned layerPattern[CSCConstants::NUM_LAYERS]) const {
0087   // 12-bit comparator code
0088   std::bitset<12> cc(code);
0089 
0090   for (unsigned ilayer = 0; ilayer < CSCConstants::NUM_LAYERS; ilayer++) {
0091     std::bitset<2> cclayer;
0092     cclayer[1] = cc[2 * ilayer + 1];
0093     cclayer[0] = cc[2 * ilayer];
0094     layerPattern[ilayer] = cclayer.to_ulong();
0095   }
0096 }
0097 
0098 //fills "code_hits" with how the code "code" would look inside the pattern
0099 int CSCPattern::recoverPatternCCCombination(
0100     const unsigned code, int code_hits[CSCConstants::NUM_LAYERS][CSCConstants::CLCT_PATTERN_WIDTH]) const {
0101   // comparator code per layer
0102   unsigned layerPattern[CSCConstants::NUM_LAYERS];
0103   getLayerPattern(code, layerPattern);
0104 
0105   // now set the bits in the pattern
0106   for (unsigned int j = 0; j < CSCConstants::NUM_LAYERS; j++) {
0107     unsigned trueCounter = 0;  //for each layer, should only have 3
0108     for (unsigned int i = 0; i < CSCConstants::CLCT_PATTERN_WIDTH; i++) {
0109       // zeros in the pattern envelope, or CC0
0110       if (!pat_[j][i]) {
0111         code_hits[j][i] = 0;
0112       }
0113       // ones in the pattern envelope
0114       else {
0115         trueCounter++;
0116         if (trueCounter == layerPattern[j])
0117           code_hits[j][i] = 1;
0118         else
0119           code_hits[j][i] = 0;
0120       }
0121     }
0122   }
0123   return 0;
0124 }
0125 
0126 int convertToLegacyPattern(const int code_hits[CSCConstants::NUM_LAYERS][CSCConstants::CLCT_PATTERN_WIDTH],
0127                            unsigned N_LAYER_REQUIREMENT);
0128 void getErrors(const vector<float>& x, const vector<float>& y, float& sigmaM, float& sigmaB);
0129 void writeHeaderPosOffsetLUT(ofstream& file);
0130 void writeHeaderSlopeLUT(ofstream& file);
0131 unsigned firmwareWord(const unsigned quality, const unsigned slope, const unsigned offset);
0132 void setDataWord(unsigned& word, const unsigned newWord, const unsigned shift, const unsigned mask);
0133 unsigned assignPosition(const float fvalue, const float fmin, const float fmax, const unsigned nbits);
0134 unsigned assignBending(const float fvalue);
0135 
0136 int CCLUTLinearFitWriter(unsigned N_LAYER_REQUIREMENT = 3) {
0137   //all the patterns we will fit
0138   std::unique_ptr<std::vector<CSCPattern>> newPatterns(new std::vector<CSCPattern>());
0139   for (unsigned ipat = 0; ipat < 5; ipat++) {
0140     newPatterns->emplace_back(ipat, CSCPatternBank::clct_pattern_run3_[ipat]);
0141   }
0142 
0143   // create output directory
0144   const std::string outdir("output_" + std::to_string(N_LAYER_REQUIREMENT) + "layers/");
0145 
0146   // output ROOT file with fit results
0147   TFile* file = TFile::Open(TString::Format("figures_%dlayers/fitresults.root", N_LAYER_REQUIREMENT), "RECREATE");
0148 
0149   TH1D* all_offsets = new TH1D("all_offsets", "All Half-strip offsets", 200, -2.5, 2.5);
0150   all_offsets->GetYaxis()->SetTitle("Entries");
0151   all_offsets->GetXaxis()->SetTitle("Half-strip offset");
0152 
0153   TH1D* all_slopes = new TH1D("all_slopes", "All slopes", 200, -3, 3);
0154   all_slopes->GetYaxis()->SetTitle("Entries");
0155   all_slopes->GetXaxis()->SetTitle("Slope [half-strips/layer]");
0156 
0157   std::map<unsigned, TH1D*> offsets;
0158   std::map<unsigned, TH1D*> slopes;
0159   std::map<unsigned, TH1D*> offsetuncs;
0160   std::map<unsigned, TH1D*> slopeuncs;
0161   std::map<unsigned, TH1D*> chi2s;
0162   std::map<unsigned, TH1D*> chi2ndfs;
0163   std::map<unsigned, TH1D*> layers;
0164   std::map<unsigned, TH1D*> legacypatterns;
0165 
0166   for (unsigned i = 0; i < 5; ++i) {
0167     std::string title("Half-strip offset, PID " + std::to_string(i));
0168     std::string name("offset_pat" + std::to_string(i));
0169 
0170     offsets[i] = new TH1D(name.c_str(), title.c_str(), 100, -2.5, 2.5);
0171     offsets[i]->GetYaxis()->SetTitle("Entries");
0172     offsets[i]->GetXaxis()->SetTitle("Half-strip offset");
0173 
0174     title = "Slope, PID " + std::to_string(i);
0175     name = "slope_pat" + std::to_string(i);
0176     slopes[i] = new TH1D(name.c_str(), title.c_str(), 100, -3, 3);
0177     slopes[i]->GetYaxis()->SetTitle("Entries");
0178     slopes[i]->GetXaxis()->SetTitle("Slope [half-strips/layer]");
0179 
0180     title = "Half-strip offset uncertainty, PID " + std::to_string(i);
0181     name = "offsetunc_pat" + std::to_string(i);
0182     offsetuncs[i] = new TH1D(name.c_str(), title.c_str(), 100, -2, 2);
0183     offsetuncs[i]->GetYaxis()->SetTitle("Entries");
0184     offsetuncs[i]->GetXaxis()->SetTitle("Half-strip offset uncertainty");
0185 
0186     title = "Slope uncertainty, PID " + std::to_string(i);
0187     name = "slopeunc_pat" + std::to_string(i);
0188     slopeuncs[i] = new TH1D(name.c_str(), title.c_str(), 100, -2, 2);
0189     slopeuncs[i]->GetYaxis()->SetTitle("Entries");
0190     slopeuncs[i]->GetXaxis()->SetTitle("Slope uncertainty [half-strips/layer]");
0191 
0192     title = "Chi2, PID " + std::to_string(i);
0193     name = "chi2_pat" + std::to_string(i);
0194     chi2s[i] = new TH1D(name.c_str(), title.c_str(), 100, 0, 100);
0195     chi2s[i]->GetYaxis()->SetTitle("Entries");
0196     chi2s[i]->GetXaxis()->SetTitle("#chi^{2}");
0197 
0198     title = "Chi2/ndf, PID " + std::to_string(i);
0199     name = "chi2ndf_pat" + std::to_string(i);
0200     chi2ndfs[i] = new TH1D(name.c_str(), title.c_str(), 40, 0, 20);
0201     chi2ndfs[i]->GetYaxis()->SetTitle("Entries");
0202     chi2ndfs[i]->GetXaxis()->SetTitle("#chi^{2}/NDF");
0203 
0204     title = "Number of layers, PID " + std::to_string(i);
0205     name = "layers_pat" + std::to_string(i);
0206     layers[i] = new TH1D(name.c_str(), title.c_str(), 7, 0, 7);
0207     layers[i]->GetYaxis()->SetTitle("Entries");
0208     layers[i]->GetXaxis()->SetTitle("Number of layers in pattern");
0209 
0210     title = "Run-1/2 pattern, PID " + std::to_string(i);
0211     name = "legacypatterns_pat" + std::to_string(i);
0212     legacypatterns[i] = new TH1D(name.c_str(), title.c_str(), 11, 0, 11);
0213     legacypatterns[i]->GetYaxis()->SetTitle("Entries");
0214     legacypatterns[i]->GetXaxis()->SetTitle("Equivalent Run-1/2 pattern number");
0215   }
0216 
0217   for (auto& p : offsets) {
0218     (p.second)->GetYaxis()->SetTitle("Entries");
0219   }
0220 
0221   //Used to calculate span of position offsets
0222   float maxOffset = -1;
0223   float maxPatt = -1;
0224   float maxCode = -1;
0225   float minOffset = -2;
0226   float minPatt = 0;
0227   float minCode = 1;
0228 
0229   // for each pattern
0230   for (auto patt = newPatterns->begin(); patt != newPatterns->end(); ++patt) {
0231     std::cout << "Processing pattern " << patt - newPatterns->begin() << std::endl;
0232 
0233     // create 3 output files per pattern: 1) position offset for CMSSW, 2) slope for CMSSW, 3) output for firmware
0234     std::cout << "Create output files for pattern " << patt->getName() << std::endl;
0235 
0236     // floating point
0237     ofstream outoffset_sw;
0238     outoffset_sw.open(outdir + "CSCComparatorCodePosOffsetLUT_pat" + patt->getName() + "_float_v1.txt");
0239     writeHeaderPosOffsetLUT(outoffset_sw);
0240 
0241     ofstream outslope_sw;
0242     outslope_sw.open(outdir + "CSCComparatorCodeSlopeLUT_pat" + patt->getName() + "_float_v1.txt");
0243     writeHeaderSlopeLUT(outslope_sw);
0244 
0245     // unsigned
0246     ofstream outoffset_sw_bin;
0247     outoffset_sw_bin.open(outdir + "CSCComparatorCodePosOffsetLUT_pat" + patt->getName() + "_v1.txt");
0248     writeHeaderPosOffsetLUT(outoffset_sw_bin);
0249 
0250     ofstream outslope_sw_bin;
0251     outslope_sw_bin.open(outdir + "CSCComparatorCodeSlopeLUT_pat" + patt->getName() + "_v1.txt");
0252     writeHeaderSlopeLUT(outslope_sw_bin);
0253 
0254     // format: [8:0] is quality, [13, 9] is bending , [17:14] is offset
0255     ofstream outfile_fw;
0256     outfile_fw.open(outdir + "rom_pat" + patt->getName() + ".mem");
0257 
0258     // pattern conversions
0259     ofstream outpatternconv;
0260     outpatternconv.open(outdir + "CSCComparatorCodePatternConversionLUT_pat" + patt->getName() + "_v1.txt");
0261     outpatternconv << "#<header> v1.0 12 32 </header>\n";
0262 
0263     // iterate through each possible comparator code
0264     for (unsigned code = 0; code < CSCConstants::NUM_COMPARATOR_CODES; code++) {
0265       if (DEBUG > 0) {
0266         cout << "Evaluating..." << endl;
0267         patt->printCode(code);
0268       }
0269       int hits[CSCConstants::NUM_LAYERS][CSCConstants::CLCT_PATTERN_WIDTH];
0270 
0271       if (patt->recoverPatternCCCombination(code, hits)) {
0272         cout << "Error: CC evaluation has failed" << endl;
0273       }
0274 
0275       //put the coordinates in the hits into a vector
0276       // x = layer, y = position in that layer
0277       vector<float> x;
0278       vector<float> y;
0279       vector<float> xe;
0280       vector<float> ye;
0281 
0282       for (int i = 0; i < CSCConstants::NUM_LAYERS; i++) {
0283         for (unsigned j = 0; j < CSCConstants::CLCT_PATTERN_WIDTH; j++) {
0284           if (hits[i][j]) {
0285             //shift to key half strip layer (layer 3)
0286             x.push_back(i - 2);
0287             y.push_back(float(j));
0288             xe.push_back(float(0));
0289             ye.push_back(HALF_STRIP_ERROR);
0290             if (DEBUG > 0)
0291               cout << "L " << x.back() << " HS " << y.back() << endl;
0292           }
0293         }
0294         if (DEBUG > 0)
0295           cout << endl;
0296       }
0297 
0298       // fit results
0299       float offset = 0;
0300       float slope = 0;
0301       float offsetunc = 0;
0302       float slopeunc = 0;
0303       float chi2 = -1;
0304       unsigned ndf;
0305       unsigned layer = 0;
0306       unsigned legacypattern = 0;
0307 
0308       // consider at least patterns with 3 hits
0309       if (x.size() >= N_LAYER_REQUIREMENT) {
0310         // for each combination fit a straight line
0311         std::unique_ptr<TGraphErrors> gr(new TGraphErrors(x.size(), &x[0], &y[0], &xe[0], &ye[0]));
0312         std::unique_ptr<TF1> fit(new TF1("fit", "pol1", -3, 4));
0313         gr->Fit("fit", "EMQ");
0314 
0315         // fit results
0316         // center at the key half-strip
0317         offset = fit->GetParameter(0) - halfpatternwidth;
0318         slope = fit->GetParameter(1);
0319         chi2 = fit->GetChisquare();
0320         ndf = fit->GetNDF();
0321         layer = ndf + 2;
0322         legacypattern = convertToLegacyPattern(hits, N_LAYER_REQUIREMENT);
0323         // mean half-strip error; slope half-strip error
0324         getErrors(x, y, offsetunc, slopeunc);
0325 
0326         // save fit results
0327         const unsigned ipat(patt - newPatterns->begin());
0328         offsets[ipat]->Fill(offset);
0329         slopes[ipat]->Fill(slope);
0330         offsetuncs[ipat]->Fill(offsetunc);
0331         slopeuncs[ipat]->Fill(slopeunc);
0332         chi2s[ipat]->Fill(chi2);
0333         chi2ndfs[ipat]->Fill(chi2 / ndf);
0334         layers[ipat]->Fill(layer);
0335         legacypatterns[ipat]->Fill(legacypattern);
0336         // all results
0337         all_offsets->Fill(offset);
0338         all_slopes->Fill(slope);
0339       }
0340 
0341       // everything is in half-strips
0342       const float fmaxOffset = 2;
0343       const float fminOffset = -1.75;
0344 
0345       // negative bending -> 0
0346       // positive bending -> 1
0347       const bool slope_sign(slope >= 0);
0348 
0349       const unsigned offset_bin = assignPosition(offset, fminOffset, fmaxOffset, 4);
0350       unsigned slope_bin = assignBending(std::abs(slope));
0351       if (slope_sign)
0352         slope_bin += 16;
0353       const unsigned fwword = firmwareWord(0, slope_bin, offset_bin);
0354 
0355       // write to output files
0356       outoffset_sw << code << " " << offset << "\n";
0357       outslope_sw << code << " " << slope << "\n";
0358       outoffset_sw_bin << code << " " << offset_bin << "\n";
0359       outslope_sw_bin << code << " " << slope_bin << "\n";
0360       outpatternconv << code << " " << legacypattern << "\n";
0361       outfile_fw << setfill('0');
0362       outfile_fw << setw(5) << std::hex << fwword << "\n";
0363 
0364       // calculate min and max codes
0365       if (layer >= N_LAYER_REQUIREMENT) {
0366         if (offset < minOffset) {
0367           minOffset = offset;
0368           minPatt = patt->getId();
0369           minCode = code;
0370         }
0371         if (offset > maxOffset) {
0372           maxOffset = offset;
0373           maxPatt = patt->getId();
0374           maxCode = code;
0375         }
0376       }
0377     }  // end loop on comparator codes
0378 
0379     // write to files
0380     outoffset_sw.close();
0381     outslope_sw.close();
0382     outoffset_sw_bin.close();
0383     outslope_sw_bin.close();
0384     outfile_fw.close();
0385   }
0386 
0387   cout << "minOffset = " << minOffset << endl;
0388   for (auto patt = newPatterns->begin(); patt != newPatterns->end(); ++patt) {
0389     if (patt->getId() == minPatt) {
0390       patt->printCode(minCode);
0391     }
0392   }
0393   cout << "maxOffset = " << maxOffset << endl;
0394   for (auto patt = newPatterns->begin(); patt != newPatterns->end(); ++patt) {
0395     if (patt->getId() == maxPatt) {
0396       patt->printCode(maxCode);
0397     }
0398   }
0399 
0400   // plot the figures
0401   for (auto& q : {offsets, slopes, offsetuncs, slopeuncs, chi2ndfs, layers, legacypatterns}) {
0402     for (auto& p : q) {
0403       TCanvas* c = new TCanvas((p.second)->GetName(), (p.second)->GetTitle(), 800, 600);
0404       c->cd();
0405       gPad->SetLogy(1);
0406       (p.second)->Draw();
0407       c->SaveAs(TString::Format("figures_%dlayers/", N_LAYER_REQUIREMENT) + TString((p.second)->GetName()) + ".pdf");
0408     }
0409   }
0410 
0411   file->Write();
0412   file->Close();
0413 
0414   return 0;
0415 }
0416 
0417 /// helpers
0418 
0419 int searchForHits(const int code_hits[CSCConstants::NUM_LAYERS][CSCConstants::CLCT_PATTERN_WIDTH],
0420                   int delta,
0421                   unsigned N_LAYER_REQUIREMENT) {
0422   unsigned returnValue = 0;
0423   unsigned maxlayers = 0;
0424   for (unsigned iPat = 2; iPat < CSCPatternBank::clct_pattern_legacy_.size(); iPat++) {
0425     unsigned nlayers = 0;
0426     for (int i = 0; i < CSCConstants::NUM_LAYERS; i++) {
0427       bool layerhit = false;
0428       for (unsigned j = 0; j < CSCConstants::CLCT_PATTERN_WIDTH; j++) {
0429         // shifted index
0430         int jdelta = j + delta;
0431         if (jdelta < 0)
0432           jdelta = 0;
0433         if (jdelta >= CSCConstants::CLCT_PATTERN_WIDTH)
0434           jdelta = CSCConstants::CLCT_PATTERN_WIDTH - 1;
0435 
0436         // do not consider invalid pattern hits
0437         if (CSCPatternBank::clct_pattern_legacy_[iPat][i][jdelta]) {
0438           // is the new pattern half-strip hit?
0439           if (code_hits[i][j]) {
0440             layerhit = true;
0441           }
0442         }
0443       }
0444       // increase the number of layers hit
0445       if (layerhit) {
0446         nlayers++;
0447       }
0448     }
0449     if (nlayers > maxlayers and nlayers >= N_LAYER_REQUIREMENT) {
0450       maxlayers = nlayers;
0451       returnValue = iPat;
0452     }
0453   }
0454   return returnValue;
0455 }
0456 
0457 // function to convert
0458 int convertToLegacyPattern(const int code_hits[CSCConstants::NUM_LAYERS][CSCConstants::CLCT_PATTERN_WIDTH],
0459                            unsigned N_LAYER_REQUIREMENT) {
0460   int returnValue = searchForHits(code_hits, 0, N_LAYER_REQUIREMENT);
0461   // try the search on a half-strip to the left
0462   if (!returnValue)
0463     returnValue = searchForHits(code_hits, -1, N_LAYER_REQUIREMENT);
0464   // try the search on a half-strip to the right
0465   if (!returnValue)
0466     returnValue = searchForHits(code_hits, 1, N_LAYER_REQUIREMENT);
0467   return returnValue;
0468 }
0469 
0470 void getErrors(const vector<float>& x, const vector<float>& y, float& sigmaM, float& sigmaB) {
0471   int N = x.size();
0472   if (!N) {
0473     sigmaM = -9;
0474     sigmaB = -9;
0475     return;
0476   }
0477 
0478   float sumx = 0;
0479   float sumx2 = 0;
0480   for (int i = 0; i < N; i++) {
0481     sumx += x[i];
0482     sumx2 += x[i] * x[i];
0483   }
0484 
0485   float delta = N * sumx2 - sumx * sumx;
0486 
0487   sigmaM = HALF_STRIP_ERROR * sqrt(1. * N / delta);
0488   sigmaB = HALF_STRIP_ERROR * sqrt(sumx2 / delta);
0489 }
0490 
0491 void writeHeaderPosOffsetLUT(ofstream& file) {
0492   file << "#CSC Position Offset LUT\n"
0493        << "#Structure is comparator code (iCC) -> Stub position offset\n"
0494        << "#iCC bits: 12\n"
0495        << "#<header> v1.0 12 32 </header>\n";
0496 }
0497 
0498 void writeHeaderSlopeLUT(ofstream& file) {
0499   file << "#CSC Slope LUT\n"
0500        << "#Structure is comparator code (iCC) -> Stub slope\n"
0501        << "#iCC bits: 12\n"
0502        << "#<header> v1.0 12 32 </header>\n";
0503 }
0504 
0505 unsigned assignPosition(const float fvalue, const float fmin, const float fmax, const unsigned nbits) {
0506   bool debug = false;
0507   unsigned value = 0;
0508   const unsigned range = pow(2, nbits);
0509   const unsigned minValue = 0;
0510   const unsigned maxValue = range - 1;
0511   const double fdelta = (fmax - fmin) / range;
0512 
0513   if (fvalue >= fmax) {
0514     value = maxValue;
0515   } else if (fvalue <= fmin) {
0516     value = minValue;
0517   } else {
0518     value = std::min(unsigned(std::ceil((fvalue - fmin) / fdelta)), maxValue);
0519   }
0520   if (debug)
0521     std::cout << "fvalue " << fvalue << " " << fmin << " " << fmax << " " << nbits << " " << value << std::endl;
0522 
0523   return value;
0524 }
0525 
0526 unsigned assignBending(const float fvalue) {
0527   bool debug = false;
0528   unsigned value = 0;
0529 
0530   // as defined in DN-19-059, section 4.8
0531   float slopes[17] = {
0532       0.0, 0.125, 0.25, 0.375, 0.5, 0.625, 0.75, 0.875, 1.0, 1.125, 1.25, 1.375, 1.5, 1.625, 1.75, 2.0, 2.5};
0533 
0534   for (unsigned i = 0; i < 16; i++) {
0535     if (fvalue >= slopes[i] and fvalue < slopes[i + 1]) {
0536       value = i;
0537     }
0538   }
0539   // overflow bin or undefined value
0540   if (fvalue >= slopes[16]) {
0541     value = 15;
0542   }
0543   return value;
0544 }
0545 
0546 unsigned firmwareWord(const unsigned quality, const unsigned slope, const unsigned offset) {
0547   /* construct fw dataword:
0548      [8:0] is quality (set all to 0 for now)
0549      [12:9] is slope value
0550      [13] is slope sign
0551      [17:14] is offset
0552   */
0553   enum Masks { OffsetMask = 0xf, SlopeMask = 0x1f, QualityMask = 0x1ff };
0554   enum Shifts { OffsetShift = 14, SlopeShift = 9, QualityShift = 0 };
0555 
0556   unsigned fwword = 0;
0557   setDataWord(fwword, quality, QualityShift, QualityMask);
0558   setDataWord(fwword, slope, SlopeShift, SlopeMask);
0559   setDataWord(fwword, offset, OffsetShift, OffsetMask);
0560 
0561   return fwword;
0562 }
0563 
0564 void setDataWord(unsigned& word, const unsigned newWord, const unsigned shift, const unsigned mask) {
0565   // clear the old value
0566   word &= ~(mask << shift);
0567 
0568   // set the new value
0569   word |= newWord << shift;
0570 }