Changeset 11917 for trunk/FACT++/gui/RawEventsViewer

Timestamp:

08/30/11 13:52:37 (13 years ago)

Author:

lyard

Message:

Added pixels 1438 and 1439

File:

• trunk/FACT++/gui/RawEventsViewer/RawEventsViewer.cc (modified) (17 diffs)

trunk/FACT++/gui/RawEventsViewer/RawEventsViewer.cc

@@ -201 +201 @@
     for (int i=0;i<ACTUAL_NUM_PIXELS;i++)
     {
-        if (i == 690 ||
-            i == 70)
-            continue;
+//        if (i == 690 ||
+//            i == 70)
+//            continue;
         glColor3fv(pixelsColor[i]);
         glLoadName(i);
@@ -217 +217 @@
     for (int i=0;i<ACTUAL_NUM_PIXELS;i++)
     {
-        if (i == 690 ||
-            i == 70)
-            continue;
+//        if (i == 690 ||
+//            i == 70)
+//            continue;
         drawHexagon(i, false);
     }
@@ -292 +292 @@
     glVertex2f(bboxMin[0] + xRange*whichSlice/(float)nRoi,
                bboxMax[1]);
+
+/*    glColor3f(0.f,0.5f,0.f);
+    for (int i=0;i<nRoi-1;i++)
+    {
+        glVertex2f(bboxMin[0] + xRange*i/(float)nRoi,
+                   bboxMin[1] + yRange*(n1mean[ i]+plusMe) /divideMe);
+        glVertex2f(bboxMin[0] + xRange*(i+1)/(float)nRoi,
+                   bboxMin[1] + yRange*(n1mean[i+1]+plusMe) /divideMe);
+    }
+*/
     glEnd();
     glEnable(GL_MULTISAMPLE);
@@ -448 +458 @@
         delete[] eventData;
         delete[] rawEventData;
+        delete[] waveLetArray;
     }
 }
@@ -461 +472 @@
         for (int j=0;j<9;j++)//for all cells of the current patch
         {
-            if (patches[i][j] == 690 ||
-                patches[i][j] == 70)
-                continue;
+//            if (patches[i][j] == 690 ||
+//                patches[i][j] == 70)
+//                continue;
             for (int k=0;k<6;k++)//for all sides of the current cell
             {
@@ -685 +696 @@
         delete[] eventData;
         delete[] rawEventData;
+        delete[] waveLetArray;
     }
     eventData = new float[(1440+160)*nRoi];
     rawEventData = new int16_t[(1440+160)*nRoi];
+    waveLetArray = new int16_t[1024*1440];
     if (!inputFile->SetPtrAddress("Data", rawEventData)){
         cout << "Missing column " << "Data" << " Aborting load..." << endl;
@@ -715 +728 @@
     rowNum = -1;
     eventStep = 1;
+
+    doMyWaveletTestPlease();
+
+
     plusEvent();
     eventStep = backupStep;
     emit newFileLoaded();
     emit signalCurrentPixel(selectedPixel);
+}
+/*
+bool RawDataViewer::doWaveLetsPlease(int givenSpan, int16_t* orig_in, int16_t* wavelets_in, bool verifyResult)
+{
+    float* interArray = new float[givenSpan];
+    float* resultArray = new float[givenSpan];
+    float* orig = new float[givenSpan];
+    float* wavelets = new float[givenSpan];
+    for (int i=0;i<givenSpan;i++)
+    {
+        if (fabs((float)(orig_in[i])) > 32767.f)
+            cout << "Input overflow: " << orig_in[i] << endl;
+        orig[i] = (float)orig_in[i];
+    }
+    for (int k=0;k<givenSpan;k++)
+        interArray[k] = orig[k];
+    int span = givenSpan/2;
+    while (span > 0)
+    {
+        for (int k=0;k<span;k++)
+        {
+            wavelets[k]        = (interArray[2*k] + interArray[2*k + 1])/2.f;
+            wavelets[k + span] = (interArray[2*k + 1] - interArray[2*k])/2.f;
+        }
+        for (int k=0;k<givenSpan;k++)
+            interArray[k] = wavelets[k];
+        span /= 2;
+    }
+    //move float results to int16_t array
+    float intScaling = 15.9f;
+    int max = 0;
+    for (int i=0;i<givenSpan;i++)
+    {
+        float cValue = intScaling*wavelets[i];//
+        if (cValue > 0)
+            cValue += 0.5f;
+        else
+            cValue -= 0.5f;
+        wavelets_in[i] = (int16_t)(cValue);
+        if (fabs(cValue) >  32767.f)
+        {
+            cout << "Overflow ! " << cValue << endl;
+        }
+        if (fabs(cValue) > fabs(max))
+            max = cValue;
+    }
+//    cout << "Max wave value: " << max << endl;
+    //result reconstruction and checking
+    if (!verifyResult)
+        return true;
+
+   for (int k=0;k<givenSpan;k++)
+   {
+       resultArray[k] = wavelets_in[k]/intScaling;
+   }
+
+   span = 1;
+   while (span < givenSpan)
+   {
+       for (int k=0;k<givenSpan;k++)
+           interArray[k] = resultArray[k];
+      for (int k=0;k<span;k++)
+       {
+          resultArray[2*k]     = (float)(((interArray[k] - interArray[k + span])*(1.f/1.f)) + 0.0f);
+          resultArray[2*k + 1] = (float)(((interArray[k] + interArray[k + span])*(1.f/1.f)) + 0.0f);
+       }
+       span *= 2;
+   }
+
+   for (int k=0;k<givenSpan;k++)
+   {
+       float plus = 0.5f;
+       if (resultArray[k] < 0)
+           plus *= -1.f;
+       if ((int)(resultArray[k]+plus) != (int)(orig_in[k]))
+       {
+           cout << "Nop, sorry: k: "  << k << " " << resultArray[k] << " " << (int)(resultArray[k]+plus) << " " << plus << " " << orig[k] << endl;
+           return false;
+       }
+   }
+   return true;
+}
+*/
+bool RawDataViewer::doWaveLetsPlease(int givenSpan, int16_t* orig_in, int16_t* wavelets_in, bool verifyResult)
+{
+    float* interArray = new float[givenSpan];
+    float* resultArray = new float[givenSpan];
+    float* orig = new float[givenSpan];
+    float* wavelets = new float[givenSpan];
+    for (int i=0;i<givenSpan;i++)
+    {
+        if (fabs((float)(orig_in[i])) > 32767.f)
+            cout << "Input overflow: " << orig_in[i] << endl;
+        orig[i] = (float)orig_in[i];
+    }
+    for (int k=0;k<givenSpan;k++)
+        interArray[k] = orig[k];
+    int span = givenSpan/2;
+    while (span > 0)
+    {
+        for (int k=0;k<span;k++)
+        {
+            wavelets[k]        = interArray[2*k];// + interArray[2*k + 1])/2.f;
+            wavelets[k + span] = interArray[2*k + 1] - interArray[2*k];//)/2.f;
+        }
+        for (int k=0;k<givenSpan;k++)
+            interArray[k] = wavelets[k];
+        span /= 2;
+    }
+    //move float results to int16_t array
+    float intScaling = 1.f;
+    int max = 0;
+    for (int i=0;i<givenSpan;i++)
+    {
+        float cValue = intScaling*wavelets[i];//
+        if (cValue > 0)
+            cValue += 0.5f;
+        else
+            cValue -= 0.5f;
+        wavelets_in[i] = (int16_t)(cValue);
+        if (fabs(cValue) >  32767.f)
+        {
+            cout << "Overflow ! " << cValue << endl;
+        }
+        if (fabs(cValue) > fabs(max))
+            max = cValue;
+    }
+//    cout << "Max wave value: " << max << endl;
+    //result reconstruction and checking
+    if (!verifyResult)
+        return true;
+
+   for (int k=0;k<givenSpan;k++)
+   {
+       resultArray[k] = wavelets_in[k]/intScaling;
+   }
+
+   span = 1;
+   while (span < givenSpan)
+   {
+       for (int k=0;k<givenSpan;k++)
+           interArray[k] = resultArray[k];
+      for (int k=0;k<span;k++)
+       {
+          resultArray[2*k]     = (float)(((interArray[k])*(1.f/1.f)) + 0.0f);
+          resultArray[2*k + 1] = (float)(((interArray[k] + interArray[k + span])*(1.f/1.f)) + 0.0f);
+       }
+       span *= 2;
+   }
+
+   for (int k=0;k<givenSpan;k++)
+   {
+       float plus = 0.5f;
+       if (resultArray[k] < 0)
+           plus *= -1.f;
+       if ((int)(resultArray[k]+plus) != (int)(orig_in[k]))
+       {
+           cout << "Nop, sorry: k: "  << k << " " << resultArray[k] << " " << (int)(resultArray[k]+plus) << " " << plus << " " << orig[k] << endl;
+           return false;
+       }
+   }
+   delete[] interArray;
+   delete[] resultArray;
+   delete[] orig;
+   delete[] wavelets;
+   return true;
+}
+
+void RawDataViewer::doWaveLetOnCurrentEventPlease()
+{
+     int16_t* origTheWayIWant = new int16_t[1024*1440];
+
+     int numPixels = 1024;
+         int leftOver = 1440-numPixels;
+         for (int k=0;k<1024;k++)
+             for (int j=0;j<1440;j++)
+             {
+                 origTheWayIWant[k*1440 + j] = rawEventData[j + k*1440];
+             }
+
+         int waveSpan = 1024;
+         int waveToPixelsRatio = 32;
+         waveSpan = waveSpan*waveToPixelsRatio;
+         int totalNumSamples = 1024*1440;
+         int j=0;
+         for (;j<totalNumSamples;j+= waveSpan)
+         {
+             if (j + waveSpan < totalNumSamples)
+             if (!doWaveLetsPlease(waveSpan, &origTheWayIWant[j], &waveLetArray[j], true))
+                 return;
+         }
+
+         while (j%1440 != 0)
+         {
+             int lastRun = 1;
+             while (lastRun < 1440 - j%1440)
+                 lastRun *= 2;
+             lastRun /= 2;
+             if (lastRun > 2)
+             {
+                 doWaveLetsPlease(lastRun, &origTheWayIWant[j], &waveLetArray[j], true);
+             }
+             else
+             {
+                 for (int l=0;l<lastRun;l++)
+                     waveLetArray[j+l] = origTheWayIWant[j+l];
+             }
+             if (!lastRun)
+                 break;
+             j += lastRun;
+         }
+
+     delete[] origTheWayIWant;
+
+}
+void RawDataViewer::doMyWaveletTestPlease()
+{
+//    cout << "Size of float: " << sizeof(float) << endl;
+    return;
+    ofstream outBin("/scratch/bin/outputBin.bin", ios_base::out | ios_base::binary);
+    ofstream outWave("/scratch/bin/outputWave.bin", ios_base::out | ios_base::binary);
+    int16_t* waveLetArray = new int16_t[1024*1440];
+    int16_t* origTheWayIWant = new int16_t[1024*1440];
+
+    for (int i=0;i<300;i++)//nRows;i++)
+    {
+        cout << '\r' <<  "Doing row " << i << " of " << nRows;
+        cout.flush();
+        inputFile->GetRow(i);
+        outBin.write((const char*)(rawEventData), 1440*1024*2);
+
+        int numPixels = 1024;
+        int leftOver = 1440-numPixels;
+        for (int k=0;k<1024;k++)
+            for (int j=0;j<1440;j++)
+            {
+                origTheWayIWant[k*1440 + j] = rawEventData[j + k*1440];
+            }
+
+        int waveSpan = 1024;
+        int waveToPixelsRatio = 32;
+        waveSpan = waveSpan*waveToPixelsRatio;
+        int totalNumSamples = 1024*1440;
+        int j=0;
+        for (;j<totalNumSamples;j+= waveSpan)//1440/waveToPixelsRatio;j++)
+        {
+            if (j + waveSpan < totalNumSamples)
+            if (!doWaveLetsPlease(waveSpan, &origTheWayIWant[j], &waveLetArray[j], true))
+                return;
+        }
+
+        while (j%1440 != 0)// < totalNumSamples)
+        {
+            int lastRun = 1;
+            while (lastRun < 1440 - j%1440)//totalNumSamples-j)
+                lastRun *= 2;
+            lastRun /= 2;
+            if (lastRun > 2)
+            {
+//                   cout << "   Doint one last run of " << lastRun << " samples" << endl;
+                doWaveLetsPlease(lastRun, &origTheWayIWant[j], &waveLetArray[j], true);
+            }
+            else
+            {
+//                    cout << " Filling in " << lastRun << " samples" << endl;
+                for (int l=0;l<lastRun;l++)
+                    waveLetArray[j+l] = origTheWayIWant[j+l];
+            }
+            if (!lastRun)
+                break;
+            j += lastRun;
+        }
+       outWave.write((const char*)(waveLetArray), 1440*1024*2);
+    }
+    outWave.close();
+    outBin.close();
+    inputFile->GetRow(0);
+
+    delete[] waveLetArray;
+    delete[] origTheWayIWant;
 }
 void RawDataViewer::openCalibFile(string& file)
@@ -851 +1148 @@
     emit signalCurrentSlice(whichSlice);
     updateGL();
+}
+
+void RawDataViewer::computePulsesStatistics()
+{
+    if (!inputFile)
+    {
+        cout << "A FITS file must be open in order to complete this operation" << endl;
+        return;
+    }
+
+
+//    for (int i=0;i<nRows;i++)//for all events
+//    {
+//        inputFile->GetRow(rowNum);
+//        for (int i=0;i<(1440+160)*nRoi;i++)
+//            eventData[i] = (float)rawEventData[i];
+
+//        for (int j=0;j<ACTUAL_NUM_PIXELS;j++)
+///        {
+    int j = selectedPixel;
+            for (int i=0;i<nRoi;i++)
+            {
+                aMeas[i] = eventData[j*1024+i];// * adcCount;
+
+            }
+            for (int i=0;i<nRoi;i++)
+            {
+                if (i==0)
+                    n1mean[i] = aMeas[i+1];
+                else
+                {
+                    if (i==1023)
+                        n1mean[i] = aMeas[i-1];
+                    else
+                        n1mean[i] = (aMeas[i-1]+aMeas[i+1])/2.f;
+                }
+            }
+            //find spike
+            for (int i=0;i<nRoi-3;i++)
+            {
+                const float fract = 0.8f;
+                float xx, xp, xpp;
+                vCorr[i] = 0;//aMeas[i];
+                xx = aMeas[i] - n1mean[i];
+                if (xx < -8.f)
+                {
+                    xp = aMeas[i+1] - n1mean[i+1];
+                    xpp = aMeas[i+2] - n1mean[i+2];
+                    if ((aMeas[i+2] - (aMeas[i] + aMeas[i+3])/2.f) > 10.f)
+                    {
+                        vCorr[i+1] = (aMeas[i] + aMeas[i+3])/2.f;
+                        vCorr[i+2] = (aMeas[i] + aMeas[i+3])/2.f;
+                        i = i+2;
+                    }
+                    else
+                    {
+                        if ((xp > -2.*xx*fract) && (xpp < -10.f))
+                        {
+                            vCorr[i+1] = n1mean[i+1];
+                            n1mean[i+2] = aMeas[i+1] - aMeas[i+3]/2.f;
+                            i++;
+                        }
+                    }
+                }
+            }
+            for (int i=0;i<nRoi;i++)
+                n1mean[i] = aMeas[i]-n1mean[i];
+ //       }
+ //   }
 }
 /************************************************************
@@ -1263 +1629 @@
  * CURRENT EVENT HAS CHANGED. ui handler
  ************************************************************/
-double xval[4096];
-double yval[4096];
+double xval[50000];
+double yval[50000];
 void UIConnector::eventChangedFromSpinner(int cEvent)
 {
@@ -1286 +1652 @@
     viewer->updateGL();
 }
-void UIConnector::currentEventHasChanged(int cEvent)
+void UIConnector::currentEventHasChanged(int )
 {
     ostringstream str;
@@ -1298 +1664 @@
         updateSpinnerDisplay = true;
     }
+
+    viewer->doWaveLetOnCurrentEventPlease();
+
         //retrieve the data that we want to display
     str.str("");
@@ -1356 +1725 @@
             pixelHistoMap[viewer->startPix[i]] = 1;
     }
+
     std::map<int, int> timeMarksMap;
     for (int i=0;i <NTMARK; i++)
@@ -1540 +1910 @@
        triggerDelayHistoItem.setSamples(xval, yval, nsamples);
 #endif
-
-//    startCellHistoZoom->setZoomBase(startCellHistoItem.boundingRect());
+       //WAVELETS HACK
+       std::map<int, int> valuesHistoMap;
+       std::map<int, int> waveletHistoMap;
+       for (int i=0;i<1024*1440;i++)
+       {
+           if (valuesHistoMap.find(viewer->rawEventData[i]) != valuesHistoMap.end())
+               valuesHistoMap[viewer->rawEventData[i]]++;
+           else
+               valuesHistoMap[viewer->rawEventData[i]] = 1;
+           if (waveletHistoMap.find(viewer->waveLetArray[i]) != waveletHistoMap.end())
+               waveletHistoMap[viewer->waveLetArray[i]]++;
+           else
+               waveletHistoMap[viewer->waveLetArray[i]] = 1;
+       }
+
+       it = valuesHistoMap.begin();
+       nsamples = 0;
+       previousValue = it->first-10;
+       cout << "Num values Original: " << valuesHistoMap.size() << endl;
+       for (unsigned int i=0;i<valuesHistoMap.size();i++)
+       {
+           if (previousValue != it->first-1)
+           {
+               xval[nsamples] = previousValue+1;
+               yval[nsamples] = 0;
+               nsamples++;
+               xval[nsamples] = it->first-1;
+               yval[nsamples] = 0;
+               nsamples++;
+           }
+           xval[nsamples] = it->first;
+           yval[nsamples] = it->second;
+           previousValue = it->first;
+           it++;
+           nsamples++;
+           xval[nsamples] = previousValue;
+           yval[nsamples] = 0;
+           nsamples++;
+           if (nsamples > 50000)
+           {
+               cout << "Error: Maximum number of samples reached for histograms. skipping what's remaining" << endl;
+               break;
+           }
+       }
+       xval[nsamples] = it==valuesHistoMap.begin() ? 0 : (--it)->first+1;
+       yval[nsamples] = 0;
+       nsamples++;
+     //  if (nsamples > 5)
+   #if QWT_VERSION < 0x060000
+          triggerDelayHistoItem.setData(xval, yval, nsamples);
+   #else
+          triggerDelayHistoItem.setSamples(xval, yval, nsamples);
+   #endif
+
+          it = waveletHistoMap.begin();
+          nsamples = 0;
+          previousValue = it->first-10;
+          cout << "Num values WaveLets: " << waveletHistoMap.size() << endl;
+          for (unsigned int i=0;i<waveletHistoMap.size();i++)
+          {
+              if (previousValue != it->first-1)
+              {
+                  xval[nsamples] = previousValue+1;
+                  yval[nsamples] = 0;
+                  nsamples++;
+                  xval[nsamples] = it->first-1;
+                  yval[nsamples] = 0;
+                  nsamples++;
+              }
+              xval[nsamples] = it->first;
+              yval[nsamples] = it->second;
+              previousValue = it->first;
+              it++;
+              nsamples++;
+              xval[nsamples] = previousValue;
+              yval[nsamples] = 0;
+              nsamples++;
+              if (nsamples > 50000)
+              {
+                  cout << "Error: Maximum number of samples reached for histograms. skipping what's remaining" << endl;
+                  break;
+              }
+          }
+          xval[nsamples] = it==waveletHistoMap.begin() ? 0 : (--it)->first+1;
+          yval[nsamples] = 0;
+          nsamples++;
+        //  if (nsamples > 5)
+      #if QWT_VERSION < 0x060000
+          startTimeMarkHistoItem.setData(xval, yval, nsamples);
+      #else
+          startTimeMarkHistoItem.setSamples(xval, yval, nsamples);
+      #endif
+
+//END OF WAVELETS HACK
+       //    startCellHistoZoom->setZoomBase(startCellHistoItem.boundingRect());
     QStack< QRectF > stack;
 //    QRectF cRectangle = boardsTimeHistoItem.boundingRect();
@@ -1645 +2108 @@
     //curve
 //    pixelValueCurveItem.setSymbol(new QwtSymbol(QwtSymbol::Cross, Qt::NoBrush, QPen(Qt::black), QSize(5,5)));
-    pixelValueCurveItem.setPen(QColor(Qt::darkGreen));
+    pixelValueCurveItem.setPen(QColor(Qt::black));
+    aMeanCurveItem.setPen(QColor(Qt::darkGreen));
+    vCorrCurveItem.setPen(QColor(Qt::red));
+    meanCurveItem.setPen(QColor(Qt::blue));
     pixelValueCurveItem.setStyle(QwtPlotCurve::Lines);
+    aMeanCurveItem.setStyle(QwtPlotCurve::Lines);
+    vCorrCurveItem.setStyle(QwtPlotCurve::Lines);
+    meanCurveItem.setStyle(QwtPlotCurve::Lines);
+
 //    pixelValueCurveItem.setCurveAttribute(QwtPlotCurve::Fitted);
     pixelValueCurveItem.attach(pixelValueCurve);
+//    aMeanCurveItem.attach(pixelValueCurve);
+    vCorrCurveItem.attach(pixelValueCurve);
+//    meanCurveItem.attach(pixelValueCurve);
 
     //FIXME delete these pointers with the destructor
@@ -1685 +2158 @@
 #endif
     }
+
+    viewer->computePulsesStatistics();
 #if QWT_VERSION < 0x060000
-       pixelValueCurveItem.setData(xval, yval, viewer->nRoi);
+    pixelValueCurveItem.setData(xval, yval, viewer->nRoi);
+    aMeanCurveItem.setData(xval, viewer->aMeas, viewer->nRoi);
+    meanCurveItem.setData(xval, viewer->n1mean, viewer->nRoi);
+    vCorrCurveItem.setData(xval, viewer->vCorr, viewer->nRoi-3);
 #else
        pixelValueCurveItem.setSamples(xval, yval, viewer->nRoi);
+       aMeanCurveItem.setSamples(xval, viewer->aMeas, viewer->nRoi);
+       meanCurveItem.setSamples(xval, viewer->n1mean, viewer->nRoi);
+       vCorrCurveItem.setSamples(xval, viewer->vCorr, viewer->nRoi-3);
 #endif
 
+
     QStack< QRectF > stack;
-    stack.push(scaleBoundingRectangle(pixelValueCurveItem.boundingRect(), 1.05f));
+    stack.push(scaleBoundingRectangle(pixelValueCurveItem.boundingRect(), 1.5f));
     curveZoom->setZoomStack(stack);
     stack.pop();
@@ -1698 +2180 @@
     hwID = pixel;
     swID = viewer->softwareMapping[pixel];
-//    cout << "here" << endl;
-/*    int channel = hwID%360;
-    channel = channel%36;
-    channel = channel%9;
-    int inter = (pixel-channel)%360;
-     inter = inter%36;
-    patchID = inter/9;
-    inter = pixel - channel - patchID*9;
-    inter = inter%360;
-    boardID = inter/36;
-    inter = pixel - channel - patchID*9 - boardID*36;
-    crateID = inter/360;
-*/
     crateID = pixel/360;
     boardID = (pixel - crateID*360)/36;
     patchID = (pixel - crateID*360 - boardID*36)/9;
-//    cout << "Value: " << hwID << " " << swID;
-//    cout << " " << crateID  << " " << boardID;
-//   cout << " " << patchID  << endl;
-//    cout << "there " << crateID << " " << boardID << " " << patchID << endl;
-//    crateID = (pixel/4)/10;
- //   boardID = (pixel/4)%10;
-//    patchID = pixel%4;
-//    cout << "there2 " << crateID << " " << boardID << " " << patchID << endl;
 
     if (HwIDBox->value() != hwID)
@@ -1741 +2202 @@
     emit updateCurrentPixelSliceValue(qstr);
 
-//    cout << "Pixel: " << viewer->softwareMapping[pixel] << " Hardware ID: " << pixel << " Patch: " << pixel%4;
-//    cout << " Board: " << (pixel/4)%10 << " Crate: " << (pixel/4)/10 << endl;
-
- //   curveZoom->setZoomBase();
-//    cout << "pixel changed ! " << pixel << endl;
 }