TES Readout Implementation

CreateCircularBuffers (NumberOfChannels,
SamplesPerWindow, TotalNumberOfWindows)
LoadFilters()
ConfigureDataChannel()

While IncomingData:
            LoadDatainBuffers()
            ApplyCrosstalkCorrection(All_Channels,CurrentWindow)

            EventMatrix = TagEvents(All_Channels)

	For each ActiveWindow : 
            WindowIncrement=+1
				Parfor each Channel:            
                	Switch (CheckEvents(EventMatrix)):
						Case No Event Before And After
						Energy = ApplyFilter(Filter00)
						Case Event Before And After
						Energy = ApplyFilter(Filter11)

						Case Event Before Only
						Energy = ApplyFilter(Filter10)

						Case Event After Only
						Energy = ApplyFilter(Filter00)


            SaveEnergyAndTimestamp(xtcFile)


		If WindowIncrement == SaveIncrement
			SaveWaveformtoXTC(xtcFile, TotalTimeWindow, AllChannels)


 


Def TagEvents:
	Parfor each channel : 
		If Slope > MinimumSlope:
			CurrentWindowEvent=True
		Else:
			CurrentWindowEvent=False
	
		EventMatrix(channel, currentWindowIndex) = CurrentWindowEvent
	Return EventMatrix

Def ApplyFilter(Filter):
	For each i, sample :
		Temp = linearFunction(sample) * sample 
		Temp = temp* Filter(i) *Filter2(i) %More filters added if necessary
		Energy = sum(temp)
		%Details on operations required for this step
		% Can we pipeline the linear function with the multiply and add?

#include <iostream>
#include <Eigen/Dense>
#include <chrono>
#include <unistd.h>
using namespace Eigen;
using namespace std;
#define FILTERS 5
#define VERBOSE 1
#define REPEAT 1000
#define SAMPLESPEREVENT 100
#define MAXWINDOW 6
#define MAXSAMPLESIZE MAXWINDOW*SAMPLESPEREVENT
#define MAXFILTERS 5
#define NUMEVENTS 1000
int evaluate_energy(const ArrayXi samples, const ArrayXXi filters, int windowsize){
    ArrayXi temp(MAXSAMPLESIZE, 1);
    int size = SAMPLESPEREVENT*windowsize;
    temp = (samples.head(size) + ArrayXi::Constant(size,52)) * 50; // linear transformation
    //cout << filters.block(0,0,size,1) << endl;
    #if FILTERS == 1
    temp = temp * filters.topRows(size);
    #elif FILTERS == 2
    temp = temp * filters.topLeftCorner(size,1) * filters.block(0,1,size,1);
    #elif FILTERS == 3
    temp = temp * filters.topLeftCorner(size,1) * filters.block(0,1,size,1) * filters.block(0,2,size,1);
    #elif FILTERS == 4
    temp = temp * filters.topLeftCorner(size,1) * filters.block(0,1,size,1) * filters.block(0,2,size,1) * filters.block(0,3,size,1);
    #elif FILTERS == 5
    temp = temp * filters.topLeftCorner(size,1) * filters.block(0,1,size,1) * filters.block(0,2,size,1) * filters.block(0,3,size,1) * filters.block(0,4,size,1);
    #endif
    int sum = temp.sum();
    //int sum = 22;
    return sum;
}
int main(int argc, char** argv)
{
    ArrayXi samples = ArrayXi::Random(MAXSAMPLESIZE);
    ArrayXXi filters = ArrayXXi::Random(MAXSAMPLESIZE, FILTERS);
    ArrayXXf tempEvents = ArrayXXf::Random(1, NUMEVENTS);
    Array<bool, 1, NUMEVENTS> isEvent;
    isEvent = tempEvents > 0.5;
    int sum;
    auto t1 = std::chrono::high_resolution_clock::now();
    for (int repeat = 0; repeat<REPEAT; repeat++){
        for (int i = 4; i<NUMEVENTS; i++){
            if (isEvent(i) == true){
                if (isEvent(i-1) == true){
                    //use filter of length 2
                    sum = evaluate_energy(samples, filters, 2);
                    #if VERBOSE > 1
                    cout << isEvent.segment<6>(i-5) << " -- " << "2 period filter used. The sum is " << sum << endl;
                    #endif
                }
                else{
                    if (isEvent(i-2) == true){
                    //use filter of length 3
                    sum = evaluate_energy(samples, filters, 3);
                    #if VERBOSE > 1
                    cout << isEvent.segment<6>(i-5) << " -- " << "3 period filter used. The sum is " << sum << endl;
                    #endif
                    }
                    else{
                        if (isEvent(i-3) == true){
                        //use filter of length 4
                        sum = evaluate_energy(samples, filters, 4);
                        #if VERBOSE > 1
                        cout << isEvent.segment<6>(i-5) << " -- " << "4 period filter used.  The sum is " << sum  << endl;
                        #endif
                        }
                        else{
                            if (isEvent(i-4) == true){
                                //use filter of length 5
                                sum = evaluate_energy(samples, filters, 5);
                                #if VERBOSE > 1
                                cout << isEvent.segment<6>(i-5) << " -- " << "5 period filter used. The sum is " << sum  << endl;
                                #endif
                            }
                            else{
                                //use single event filter
                                sum = evaluate_energy(samples, filters, 1);
                                #if VERBOSE > 1
                                cout << isEvent.segment<6>(i-5) << " -- " << "Single period filter used. The sum is " << sum  << endl;
                                #endif
                            }
                        }
                    }
                }
            }
        }
    }

    auto t2 = std::chrono::high_resolution_clock::now();
    std::chrono::duration<double, std::milli> fp_ms = t2 - t1;
    #if VERBOSE > 0
    cout << "The sum is : " << sum << endl;
    cout << "The process took " << fp_ms.count() << "ms" << endl;
    #endif
  //Multiply vectors then sum the coefficients for inner product.
}