Changeset 12471 for trunk/Mars


Ignore:
Timestamp:
11/09/11 16:44:53 (13 years ago)
Author:
tbretz
Message:
Changed the Apply function to support the time marker channels; added the time marker channels to reading the FITS file; removed some old obsolete comments; adapted AddAbs and others to the usage of any other ROI than 1024
File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/Mars/mcore/DrsCalib.h

    r12423 r12471  
    1818    std::vector<int64_t> fSum2;
    1919
    20     /*
    21     vector<map<int16_t, uint8_t> > fMedian;
    22     vector<int16_t> fResult;
    23      */
    2420public:
    2521    DrsCalibrate() : fNumEntries(0), fNumSamples(0), fNumChannels(0) { }
    … …  
    4137        fSum.resize(samples*channels);
    4238        fSum2.resize(samples*channels);
    43 
    44         //fMedian.resize(40*9*1024);
    45     }
    46 /*
    47                 continue;
    48 
    49                 if (ch<40*9)
    50                 {
    51                     fMedian[abs][val[rel]]++;
    52 
    53                     int n= 8;
    54                     if (fNumEntries>0 && fNumEntries%(n*100)==0)
    55                     {
    56                         fResult.resize(40*9*1024);
    57 
    58                         uint16_t sum = 0;
    59                         for (map<int16_t, uint8_t>::const_iterator it=fMedian[abs].begin();
    60                              it!=fMedian[abs].end(); it++)
    61                         {
    62                             map<int16_t, uint8_t>::const_iterator is = it;
    63                             is++;
    64                             if (is==fMedian[abs].end())
    65                                 break;
    66 
    67                             sum += it->second;
    68 
    69                             double med = 0;
    70 
    71                             med += int64_t(it->second)*int64_t(it->first);
    72                             med += int64_t(is->second)*int64_t(is->first);
    73                             med /= int64_t(it->second)+int64_t(is->second);
    74 
    75                             if (sum==n*50)
    76                             {
    77                                 fResult[abs] = it->first;
    78                                 cout << ch << " MED+(50): " << it->first << endl;
    79                             }
    80                             if (sum<n*50 && sum+is->second>n*50)
    81                             {
    82                                 fResult[abs] = med;
    83                                 cout << ch << " MED-(50): " << med << endl;
    84                             }
    85                         }
    86                         cout << ch << " AVG=" << float(fSum[abs])/fNumEntries << endl;
    87                         fMedian[abs].clear();
    88                     }
    89                 } // -2029 -2012 -2003 -1996 -1990 // -1955
    90                 */
     39    }
    9140
    9241    void AddRel(const int16_t *val, const int16_t *start)
    … …  
    9847                continue;
    9948
    100             const size_t pos = ch*fNumSamples;
    101 
    102             for (size_t i=0; i<fNumSamples; i++)
     49            const size_t pos = ch*1024;
     50
     51            for (size_t i=0; i<1024; i++)
    10352            {
    10453                // Value is relative to trigger
    10554                // Abs is corresponding index relative to DRS pipeline
    10655                const size_t rel = pos +  i;
    107                 const size_t abs = pos + (spos+i)%fNumSamples;
     56                const size_t abs = pos + (spos+i)%1024;
    10857
    10958                const int64_t v = val[rel];
    … …  
    14594        fNumEntries++;
    14695    }
    147 
     96    /*
    14897    void AddAbs(const int16_t *val,    const  int16_t *start,
    14998                const int32_t *offset, const uint32_t scale)
    … …  
    173122
    174123        fNumEntries++;
     124    }*/
     125
     126    void AddAbs(const int16_t *val,    const  int16_t *start,
     127                const int32_t *offset, const uint32_t scale)
     128    {
     129        // 1440 without tm, 1600 with tm
     130        for (size_t ch=0; ch<fNumChannels; ch++)
     131        {
     132            const int16_t spos = start[ch];
     133            if (spos<0)
     134                continue;
     135
     136            const size_t pos = ch*fNumSamples;
     137            const size_t drs = ch>1439 ? ((ch-1440)*9+8)*1024 : ch*1024;
     138
     139            for (size_t i=0; i<fNumSamples; i++)
     140            {
     141                // Value is relative to trigger
     142                // Offset is relative to DRS pipeline
     143                // Abs is corresponding index relative to DRS pipeline
     144                const size_t rel = pos +  i;
     145                const size_t abs = drs + (spos+i)%1024;
     146
     147                const int64_t v = int64_t(val[rel])*scale-offset[abs];
     148
     149                fSum[rel]  += v;
     150                fSum2[rel] += v*v;
     151            }
     152        }
     153
     154        fNumEntries++;
     155    }
     156
     157    static void ApplyCh(float *vec, const int16_t *val, int16_t start, uint32_t roi,
     158                        const int32_t *offset, const uint32_t scaleabs,
     159                        const int64_t *gain,   const uint64_t scalegain)
     160    {
     161        if (start<0)
     162        {
     163            memset(vec, 0, roi);
     164            return;
     165        }
     166
     167        for (size_t i=0; i<roi; i++)
     168        {
     169            // Value is relative to trigger
     170            // Offset is relative to DRS pipeline
     171            // Abs is corresponding index relative to DRS pipeline
     172            const size_t abs = (start+i)%1024;
     173
     174            const int64_t v =
     175                + int64_t(val[i])*scaleabs-offset[abs]
     176                ;
     177
     178            const int64_t div = gain[abs];
     179            vec[i] = double(v)*scalegain/div;
     180        }
    175181    }
    176182
    … …  
    191197            // Offset is relative to DRS pipeline
    192198            // Abs is corresponding index relative to DRS pipeline
    193             const size_t rel = i;
    194199            const size_t abs = (start+i)%1024;
    195200
    196201            const int64_t v =
    197                 + (int64_t(val[rel])*scaleabs-offset[abs])*scalerel
    198                 - trgoff[rel]
     202                + (int64_t(val[i])*scaleabs-offset[abs])*scalerel
     203                - trgoff[i]
    199204                ;
    200205
    201206            const int64_t div = gain[abs]*scalerel;
    202             vec[rel] = double(v)*scalegain/div;
    203         }
    204     }
    205 
    206     static void Apply(float *vec, const int16_t *val, const int16_t *start, uint32_t roi,
    207                       const int32_t *offset, const uint32_t scaleabs,
    208                       const int64_t *gain,   const uint64_t scalegain,
    209                       const int64_t *trgoff, const uint64_t scalerel)
    210     {
    211         for (size_t ch=0; ch<1440; ch++)
    212         {
    213             const size_t pos = ch*roi;
    214             const size_t drs = ch*1024;
    215             ApplyCh(vec+pos, val+pos, start[ch], roi,
    216                     offset+drs, scaleabs,
    217                     gain  +drs, scalegain,
    218                     trgoff+drs, scalerel);
     207            vec[i] = double(v)*scalegain/div;
    219208        }
    220209    }
    … …  
    282271    void GetSampleStats(float *ptr, float scale) const
    283272    {
     273        const size_t sz = fNumSamples*fNumChannels;
     274
    284275        if (fNumEntries==0)
    285276        {
    286             memset(ptr, 0, sizeof(float)*1024*1440*2);
     277            memset(ptr, 0, sizeof(float)*sz*2);
    287278            return;
    288279        }
    … …  
    293284        while (it!=fSum.end())
    294285        {
    295             *ptr             = scale*double(*it)/fNumEntries;
    296             *(ptr+1024*1440) = scale*sqrt(double(*i2*int64_t(fNumEntries) - *it * *it))/fNumEntries;
     286            *ptr      = scale*double(*it)/fNumEntries;
     287            *(ptr+sz) = scale*sqrt(double(*i2*int64_t(fNumEntries) - *it * *it))/fNumEntries;
    297288
    298289            ptr++;
    … …  
    369360
    370361    uint32_t fStep;
     362    uint16_t fRoi;   // Region of interest for trgoff
     363    uint16_t fNumTm; // Number of time marker channels in trgoff
     364
     365//    uint16_t fDAC[8];
    371366
    372367    DrsCalibration() :
    373         fOffset(1440*1024, 0),
    374         fGain(1440*1024, 4096),
    375         fTrgOff(1440*1024, 0),
     368        fOffset  (1440*1024, 0),
     369        fGain    (1440*1024, 4096),
     370        fTrgOff  (1600*1024, 0),
    376371        fNumOffset(1),
    377372        fNumGain(2000),
    … …  
    387382        fOffset.assign(1440*1024, 0);
    388383        fGain.assign  (1440*1024, 4096);
    389         fTrgOff.assign(1440*1024, 0);
     384        fTrgOff.assign(1600*1024, 0);
    390385
    391386        fNumOffset = 1;
    … …  
    424419            std::ostringstream msg;
    425420            msg << "Reading " << str << " failed: Number of rows in table is not 1.";
    426             return msg.str();
    427         }
    428 
    429         vec.resize(1440*1024*6+3);
    430 
    431         float *base = vec.data();
    432 
    433         file.SetPtrAddress("RunNumberBaseline",      base,   1);
    434         file.SetPtrAddress("RunNumberGain",          base+1, 1);
    435         file.SetPtrAddress("RunNumberTriggerOffset", base+2, 1);
    436         file.SetPtrAddress("BaselineMean",           base+0*1024*1440+3, 1024*1440);
    437         file.SetPtrAddress("BaselineRms",            base+1*1024*1440+3, 1024*1440);
    438         file.SetPtrAddress("GainMean",               base+2*1024*1440+3, 1024*1440);
    439         file.SetPtrAddress("GainRms",                base+3*1024*1440+3, 1024*1440);
    440         file.SetPtrAddress("TriggerOffsetMean",      base+4*1024*1440+3, 1024*1440);
    441         file.SetPtrAddress("TriggerOffsetRms",       base+5*1024*1440+3, 1024*1440);
    442 
    443         if (!file.GetNextRow())
    444         {
    445             std::ostringstream msg;
    446             msg << "Reading data from " << str << " failed.";
    447421            return msg.str();
    448422        }
    … …  
    452426        fNumGain   = file.GetUInt("NBGAIN");
    453427        fNumTrgOff = file.GetUInt("NBTRGOFF");
    454 
     428        fRoi       = file.GetUInt("NROI");
     429        fNumTm     = file.GetUInt("NTM");
     430/*
     431        fDAC[0]    = file.GetUInt("DAC_A");
     432        fDAC[1]    = file.GetUInt("DAC_B");
     433        fDAC[4]    = file.GetUInt("DAC_C");
     434*/
     435        vec.resize(1440*1024*4 + (1440+fNumTm)*fRoi*2 + 4);
     436
     437        float *base = vec.data();
     438
     439        file.SetPtrAddress("RunNumberBaseline",      base,   1);
     440        file.SetPtrAddress("RunNumberGain",          base+1, 1);
     441        file.SetPtrAddress("RunNumberTriggerOffset", base+2, 1);
     442        file.SetPtrAddress("BaselineMean",           base+4+0*1024*1440,                           1024*1440);
     443        file.SetPtrAddress("BaselineRms",            base+4+1*1024*1440,                           1024*1440);
     444        file.SetPtrAddress("GainMean",               base+4+2*1024*1440,                           1024*1440);
     445        file.SetPtrAddress("GainRms",                base+4+3*1024*1440,                           1024*1440);
     446        file.SetPtrAddress("TriggerOffsetMean",      base+4+4*1024*1440,                           fRoi*1440);
     447        file.SetPtrAddress("TriggerOffsetRms",       base+4+4*1024*1440+   fRoi*1440,              fRoi*1440);
     448        file.SetPtrAddress("TriggerOffsetTMMean",    base+4+4*1024*1440+ 2*fRoi*1440,              fRoi*fNumTm);
     449        file.SetPtrAddress("TriggerOffsetTMRms",     base+4+4*1024*1440+ 2*fRoi*1440+ fRoi*fNumTm, fRoi*fNumTm);
     450
     451        if (!file.GetNextRow())
     452        {
     453            std::ostringstream msg;
     454            msg << "Reading data from " << str << " failed.";
     455            return msg.str();
     456        }
     457/*
     458        fDAC[2] = fDAC[1];
     459        fDAC[4] = fDAC[1];
     460
     461        fDAC[5] = fDAC[4];
     462        fDAC[6] = fDAC[4];
     463        fDAC[7] = fDAC[4];
     464*/
    455465        fOffset.resize(1024*1440);
    456466        fGain.resize(1024*1440);
    457         fTrgOff.resize(1024*1440);
     467
     468        fTrgOff.resize(fRoi*(1440+fNumTm));
    458469
    459470        // Convert back to ADC counts: 256/125 = 4096/2000
    … …  
    461472        for (int i=0; i<1024*1440; i++)
    462473        {
    463             fOffset[i] = fNumOffset           *256*base[i+1024*1440*0+3]/125;
    464             fGain[i]   = fNumOffset*fNumGain  *256*base[i+1024*1440*2+3]/125;
    465             fTrgOff[i] = fNumOffset*fNumTrgOff*256*base[i+1024*1440*4+3]/125;
    466         }
     474            fOffset[i] = fNumOffset         *256*base[i+1024*1440*0+4]/125;
     475            fGain[i]   = fNumOffset*fNumGain*256*base[i+1024*1440*2+4]/125;
     476        }
     477
     478        for (int i=0; i<fRoi*1440; i++)
     479            fTrgOff[i] = fNumOffset*fNumTrgOff*256*base[i+1024*1440*4+4]/125;
     480
     481        for (int i=0; i<fRoi*fNumTm; i++)
     482            fTrgOff[i+1440*fRoi] = fNumOffset*fNumTrgOff*256*base[i+1024*1440*4+2*fRoi*1440+4]/125;
     483
    467484
    468485        // DAC:  0..2.5V == 0..65535
    … …  
    500517    bool IsValid() { return fStep>2; }
    501518
    502     void Apply(float *vec, int16_t *val, const int16_t *start, uint32_t roi)
    503     {
    504         DrsCalibrate::Apply(vec, val, start, roi,
    505                             fOffset.data(), fNumOffset,
    506                             fGain.data(),   fNumGain,
    507                             fTrgOff.data(), fNumTrgOff);
    508 
     519    bool Apply(float *vec, const int16_t *val, const int16_t *start, uint32_t roi)
     520    {
     521        if (roi!=fRoi)
     522        {
     523            for (size_t ch=0; ch<1440; ch++)
     524            {
     525                const size_t pos = ch*roi;
     526                const size_t drs = ch*1024;
     527
     528                DrsCalibrate::ApplyCh(vec+pos, val+pos, start[ch], roi,
     529                                      fOffset.data()+drs, fNumOffset,
     530                                      fGain.data()  +drs, fNumGain);
     531            }
     532
     533            return false;
     534        }
     535
     536        for (size_t ch=0; ch<1440; ch++)
     537        {
     538            const size_t pos = ch*fRoi;
     539            const size_t drs = ch*1024;
     540
     541            DrsCalibrate::ApplyCh(vec+pos, val+pos, start[ch], roi,
     542                                  fOffset.data()+drs, fNumOffset,
     543                                  fGain.data()  +drs, fNumGain,
     544                                  fTrgOff.data()+pos, fNumTrgOff);
     545        }
     546
     547        for (size_t ch=0; ch<fNumTm; ch++)
     548        {
     549            const size_t pos = (ch+1440)*fRoi;
     550            const size_t drs = (ch*9+8)*1024;
     551
     552            DrsCalibrate::ApplyCh(vec+pos, val+pos, start[ch], roi,
     553                                  fOffset.data()+drs, fNumOffset,
     554                                  fGain.data()  +drs, fNumGain,
     555                                  fTrgOff.data()+pos, fNumTrgOff);
     556        }
     557
     558        return true;
    509559    }
    510560};
Note: See TracChangeset for help on using the changeset viewer.