M@TE Project

{{>toc}}

# Scheduling

(FACT++: makeschedule)

The source scheduling is done mostly automatically, based on some simple rules which try to optimize the total sensitivity of the telescope. The idea behind is to optimize the schedule according to the visibility of sources and create an unbiased scheduled which does not depend on a manual scheduler.

The schedule is always calculated three days in advance depending on whether a shift is scheduled (a shifter is assigned in the shift calendar)

The shift calendar is found here:

https://www.fact-project.org/shift/

An editor to updaye the schedule manually, e.g. in case of ToO observations, is available here:

https://www.fact-project.org/schedule/

## Basic schedule

The schedule is based on a user defined list of sources. A global limit and a source specific limit can be defined on the maximum allowed zenith angle and the maximum predicted current. For each point during the night which fulfil these requirements, a relative threshold is predicted. The relative threshold $T$ is calculated as

{{{#!Latex

T=f\cdot\cos(\theta)^{-2.664}\cdot\left(\frac{I}{6.2\mu\text{A}}\right)^{0.394}

}}}

(THIS NEEDS REVIEW... WE KNOW BETTER ALREADY)

Here, $f$ is an arbitrary scale factor which is meant to increase or decrease the priority of a source in the schedule and to account for their different average flux, $I$ is the predicted current and $\theta$ the zenith angle.

The predicted current $I$ is calculated as

(FACT::Prediction)

{{{#!Latex

\frac{I}{\mu\text{A}} = 6.2\cdot 95.7\M^{2.63}\cdot \sin\theta_\text{moon}^{0.6}\cdot \left(\frac{d_\text{ME}}{384400\text{km}}\right)^{-2}\cdot \text{e}^{0.67\cdot \cos\alpha^4} + \text{e}^{-97.8+105.8\cdot\sin\beta^2} 

}}}

Here, $M$ is the illuminated fraction of the moon disc, $\theta_\text{moon}$ is the altitude of the moon, while altitudes below horizon are considered 0°. $d_\text{ME}$ is the distance between moon and Earth. The angle $\alpha$ defined the angular separation of the observed source from the moon and the angle $\beta$ is altitude of the Sun (which is negative as the Sun is supposed to be below horizon during observations).

The predictions are also plotted here: https://fact-project.org/dch/scheduling.php (add REFs)

From all sources, over the whole night, always the one with the lowest relative threshold $T$ is selected to be observed. As this can result in some unwanted artifacts like sources observed for only a few second, this schedule is revided in further steps.

## Beginning and end of the night

If the first source of the night is found to have an observation time less than 40min and it can be extended to 40min without violating any of the observation limits, it is extended to 40min. If this is not possible, it is checked if the second scheduled source can be extended to the beginning of the night without violating any observation limit. If this is not possible, then the scheduler tries to remove the second source from the schedule and extend the first source accordingly. If all that fails, the scheduled is kept as is.

An identical procedure is applied at the end of the night, where ''first'' mist be replaced by ''last'' and ''second'' by ''second last''.

## Intermediate Source

The so obtained schedule might still show unwanted artifacts.

Therefore, for the first source scheduled during the night, which is scheduled with less than 40min, is processed. If the source is scheduled in between two no-observations (SLEEP), no action is taken. If the scheduled source neither comes after a no-observation (Example: NOOBS...SOURCE...OTHER), nor before a no-observation (OTHER...SOURCE...NOOBS), nor it is included between observations of an identical source (SAME...SOURCE...SAME), then the algorithm determines in steps of one minute, the point in time after the start of its observation at which the relative threshold of the following source falls below the relative threshold of the previous source. If extending the surrounding sources to that point and removing the intermediate source from the schedule does violate any observation limit, the schedule is kept. Otherwise, the intermediate source is removed and the surrounded sources extended accordingly. This process is iteratively repeated until no further rescheduling is possible anymore.

## Mini Sources

As this procedue can still lead to very short observations which are not reasonable and not wanted, all scheduled sourced which are scheduled now with an observation time of less than 5min are removed from the schedule and replaced by a SLEEP during which the telescope does not record any data.

## Result

The current setup is:

~~~

Global maximum current: 90 uA/pix

Global zenith distance: 75 deg 

Nsources = 12 

[Mrk 421] 

[Mrk 501] 

[1ES 2344+51.4] 

[PKS 2155-304] Penalty=2.5 

[Crab] Zd<30 

[H 1426+428] Penalty=2.5 

[1ES 1959+650] Penalty=0.6 

[1ES 1218+304] Penalty=2.5 

[IC 310] Penalty=2.5 

[1H0323+342] 

[PKS 0736+01] Penalty=2.5 

[V404 Cyg] Penalty=2.5 

~~~

The schedule produced by the auto-scheduler is available at https://www.fact-project.org/dch/obs_summary.php. In addition the web-page also shows the observation eventually carried out in comparison.

# DRS Calibration

(MDrsCalibApply: MRawEvtData -> MPedestalSubtractedEvt)

## Amplitude Calibration

1. An offset calculated for each physical cell is subtracted. The offset is obtained as the average calculated from 1000 pedestal events.

2. The contents of each baseline-subtracted physical cell is multiplied by its individual gain. The gain is calculated as the average response of each physical cell to an applied constant voltage of 896.47mV. For this 1000 events are recorded.

3. An offset calculated for each logical cell is subtracted. The offset is obtained as the average calculated from 1000 pedestal events.

## Jumps correction

(MDrsCalibrate::!CorrectStep)

The baseline changes at the edged of each readout-window. That means that if consecutive readouts happen with overlapping readout windows, a sudden jump of the baseline in the readout window is visible. If the readout window is from cell i to i+r where r is the region of interest, the jumps of the baseline (only visible during the following readout) happen between the samples i+2 and i+r+9 and their following sample.

The amplitude of this effect decreases exponentially with time. Therefore, the jump correction is only applied for the past five events. For an average rate below the readout limit of the DAQ system (~200Hz) this should leave only jumps with an amplitude below 1 ADC count. Note that for some cases, the jumps might lay totally outside of the current readout window and thus their height can not be determined easily. In this case they are ignored. 

The jumps induced from all past five events are corrected individually, with the smallest (oldest) first. Upward and downward jumps (beginnin and end of readout window) are treated seperately. As the amplitude of the jump is the same in all 9 channels of each DRS chip, an average jump height <h_ch> is calculated for all steps within the range of the current readout window from the first eight channels ch=(0..7) of each chip. As the last channel is known to be slightly different in its behavior, it is excluded. To avoid a bias in the calculation from spikes or fast rising signals, the average jump height obtained from all 160 DRS4 chips is averaged H=<<h_ch>_drs>. From carefully studies, it was obtained that if their standard deviation of H exceeds 5 there is a bias in the data. The studies also showed that of 10% of the chips are excluded, a reasonable jump height can be determined. 

In a last step, the obtained average jump height is applied to all channels. While for positive (upward) jumps, all following samples are corrected, for negative (downward) jumps all previous samples are corrected.

## Spike removal

(DrsCalibrate::RemoveSpikes3, originally developped by Werner Lustermann)

If the difference between a sample i and the average of its two neighbors (i-1 and i+1) is smaller than -10 and for the following sample (i+1), this difference has opposite sign and is at least 60% larger in absolute value, it is considered a single sample wide spike and the sample i is replaced by the average of its two neighbors. 

If this difference is smaller than -5 for sample i, and it was not detected as a single spike, it is checked for a 2-sample wide spike. If this differences is larger than 5 for its two following samples i+1 and i+2, the samples i+1 and i+2 are considered a double-spike and replaced by 1/3*[i+3] + 2/3*[i] and 2/3*[i+3] + 1/3*[i], respectively.

## Sliding average

Not applied. Turned off by default.

# Saturation Treatment

(MTreatSaturation: MPedestalSubtractedEvt -> MPedestalSubtractedEvt)

Saturation is checked between sample 50 and 225 and is assumed if the maximum sample is this range is found to exceed 1900 (=1.9V). Using a 3rd-order spline, the crossing of 1.8V before and after the saturation point is determined. If the distance between both points is less than three samples, no action is taken. By a phenomenological approach (trained in ideal pulses from our paper), the amplitude and arrival time of the original pulse is estimated. Samples starting from 5 samples before the leading 1.8V crossing to 5 samples after the falling 1.8V crossing are then replaced with an ideal pulse calculated with the pulse shape from out paper. To get slightly more accurate result, a baseline is taken into account which is calculated as the average from the samples 10 to 50.

# Filter Data

(MFilterData: MPedestalSubtractedEvt -> MPedestalSubtractedEvt)

Each waveform is filtered with a 14-sample long, weighted sliding average with the following weights.

[[Image(Filter.png)]]

Weights re-sample the integral of the leading edge of our pulse (To be checked!) and are normalized such that the integral is 1. Note that this treatment logically turns the maximum determined during the signal extraction into the integral over the leading edge (is that correct?).

The filter is applied to the whole waveform.

# Software Trigger

(MSoftwareTriggerCalc: MPedestalSubtractedEvt -> MSoftwareTrigger)

The purpose of the software trigger is to mimic the real trigger with an arbitrary threshold which is high enough to keep a sample of pure shower events. The software trigger sums the waveforms of all nine pixels in a trigger patch (faulty pixels are excluded -- maybe this is not up to date!). From each summed sample, half the amplitude of the 15 sample previous to the sample is subtracted. This corresponds to the electronic clamping with a delay line in the trigger electronics. From these values, the global average and the global maximum of the minimum of each four consecutive values is obtained. The average, the global maximum and the position of the maximum is stored.

# Signal Extraction

(MExtractFact: MPedestalSubtractedEvt -> MExtractedSignalCam)

Signal extraction is applied between sample 25 and 225. The signal is extracted calculating a 2nd-order polynomial through the logarithm of each three consecutive samples. This is logically the same as fitting a Gaussian to each three points but uses an analytical approach. From all each three-points which have the global maximum of the polynomial in between, the one with the highest maximum is selected. This global maximum is defined as the signal of the pixel and the pulse position. The arrival time is defined searching backward from the pulse position for a maximum of 15 sample (7.5 ns). While every two consecutive samples are linearly interpolated, the point is searched where the leading samples cross half the height of the maximum for the first time. This point is called the arrival time of the pulse. If no leading edge is found, a uniform random number in the extraction window is assigned as arrival time. In addition, the difference between the position of the maximum and the half-height leading edge is stored.

# Amplitude Calibration

(MCalibrateFact: MExtractedSignalCam -> MSignalCam)

All pixel signals are multiplied with a factor 0.1024, to convert them to photon-equivalent units.

# Timing Calibration

(MCalibrateDrsTimes: MArrivalTimeCam -> MSignalCam)

The position of the pulse is corrected for the unequal sampling times of each sample introduced by the DRS4. The calibrated arrival time is the difference between the absolute (calibrated) time of the cell corresponding to the signal arrival time and the start cell. An additional delay obtained from the analysis of Muon rings is applied. The times are converted to nano-seconds.

# Bad Pixel Treatment

(MBadPixelsTreat: MSignalCam -> MSignalCam)

For pixels which are flagged as bad, their signal is calculated as the average signal of their neighboring pixels, if at least three valid pixels exist. The minimum and maximum arrival time of all valid neighbor pixels is determined, as well as the number and average arrival times of all pairs of neighbors in a close package configuration which have an arrival time difference of less than 3ns. If less than three neighbors belong to such pairs, the pixels is considered not to be part of a shower and the arrival time is replaces by a uniform random number between the obtained min and max arrival time. If more than three neighbors belong to such a group, the pixel is considered part of a shower and its arrival time is replaced by the obtained average. Pixels which could not be treated are flagged as unusable and are ignored in the following.

# Image Cleaning

(MImageCleanTime: MSignalCam -> MSignalCam)

In a first step islands are defined as groups of pixels in which each neighbor has an arrival time difference to at least one of its neighbors of less than $\Delta T/\Delta d=18.75ns/deg$. Note that also every pixel which does not have any neighbor within this time range forms a single-pixel sized island. For each of those individual islands, the earliest and latest arrival time within the island is obtained. If this range overlaps with the range of a neighboring island, both islands are fused into one larger island. This step is repeated until no further fusion is possible anymore. The size of an island is then defined as the sum of the signals of all pixels in the island. All pixels of an island survive image cleaning if the number of pixels N in the island and its size S exceed a user-defined threshold. Currently $N>0$ and $S>25$. The number of islands, size of the largest island and the residual size is kept. (How to handle negative sizes?)  

The remaining pixels which survived the cleaning step are called the image.

(MC Comparison: Emc, !NumPixels, How to prove that the increase of pixels is not just noise?)

# Image Parameters

## Classical "Hillas" Parameters

(MHillas::Calc: MSignalCam -> MHillas)

The basic image parameters are just a simple statistical analysis of the light distribution in the camera. The image parameters '''MeanX/MeanY''' [mm] and '''!Width/Length''' [mm] are calculated from all pixels which survived cleaning, are not marked to be ignored and have a positive signal as the first and second moment of the signal distribution in the camera (classically known as Mean and Standard Deviation). The axis along which the second moment Length is calculated is defined such that Length is maximized and hence Width minimized. The corresponding orientation of the Length-axis in the camera is called the angle '''Delta''' [rad] and referred to as ''image axis'' in the following. Delta is defined w.r.t. the x-axis and defined between -PI/2 and PI/2. As a by-product and to increase precision of further calculations, sine and cosine of the angle Delta is stored as well.

The summed signal in those pixels is called '''Size''' [pe]. Images in which the total Size is 0 or in which only one or two pixels fulfill the mentioned criteria are discarded. This parameters are also known as ''Hillas Parameters''. Very often these statistical parameters are represented by an ellipse drawn around MeanX/MeanY oriented with the angle Delta and the major any minor axis Length and Width. The area spanned by an ellipse with the major axis Width and Length is called '''Area''' [mm^2].

## Extended "Hillas" Parameters

(MHillasExt::Calc: MSignalCam -> MHillasExt)

As an extension to the classical Hillas Parameters, the third moments along the image axis, and perpendicular to it, are calculated as '''M3Long''' [mm] and '''M3Trans''' [mm]. Positive values refer to a shift in the direction of the positive x-axis.

## Timing Parameters

To account for the time development within the shower, all arrival times in the image are projected along the image axis and perpendicular to it and a linear regression on this data is calculated. The resulting slopes along the shower axis and perpendicular are called '''!SlopeLong''' [ns/mm] and '''!SlopeTrans''' [ns/mm]. Positive values refer to a shift in the direction of the positive x-axis.

(What might be interesting to test is if the separation between the position of the mean arrival time and the position of the mean signal gives separation power)

## Leakage Parameter

(MNewImagePar::Calc: MSignalCam -> MNewImagePar)

The '''Leakage''' parameter "Leakage1" describes the fraction of the signal which is located in edge pixels of the camera. Edge pixels are all pixels which have less than six direct neighbors. The number of neighbors does not take into account whether the are marked as broken or not.

## Unit conversion

To convert from camera coordinates in milli-meters to angular distances, the small-angle approximation is used. The conversion fact for a focal distance of 4889mm is ~0.0117°/mm (0.0117193246260285378).

# Source Position

## Calculation

(MSrcPosCalc: MPointingPos, MSourcePos -> MSrcPosCam)

The source position in the camera is calculated as the reflection of a single ray in the center of the reflector (the so called master-ray) on a planar surface at distance of 4889mm. For the calculation, perfect tracking is assumed and the source position is calculated from the right-ascension and declination of the observed source and the pointing position (rad/dec) at the time of each event.

## Abberation Correction

(MHillasSrc)

Due to abberation the reflected image is distorted. With Monte Carlo simulations using the existing measures reflector setup, it has been verified that image abberation shifts the center of the reflected light distribution of a point source 2% outwards. This is equal to a scale of the coordinate system. To account for this shift, the calculated source position is shifted outwards by 2% as well.

## Parameter Sign

To calculate the sign of M3Long/M3Trans and !SlopeLong/SlopeTrans w.r.t. to the source position, the source position is projected onto the image axis and perpendicular to it, respectively, with the center-of-gravity of the light distribution (MeanX/MeanY) as coordinate origin. The correct sign is then the original sign times the sign of the projection. The resulting sign then denotes whether the third moment and the time slope are pointing towards (negative) or away (positive) from the source position in the camera.

# Disp Calculation

(MFMagicCuts)

The source position is estimated using the so called Disp method (reference) in which the distance of the source from the center-of-gravity along the major axis is estimated from the geometrical distortion of the image. This distortion originates from viewing the shower under a certain inclination angle which results in an elliptic image. The estimated distance is called Disp.

Classically, the Disp-paramater P is calculated as

{{{#!Latex

P = \xi\cdot(1-\frac{W}{L})

}}}

For gammas from the source, P [°] should be identical to the distance D [here in deg] between source position and center-of-gravity. This can be used to train the method and to obtain xi. As this correlation is only relevant for gammas, ideally the obtained distance is corrected for misalignment to fall inside the PSF. Therefore, xi is optimized to obtain a minimal \Delta\xi

{{{#!Latex

\Delta\xi = \xi-\frac{D\cdot\sec\alpha}{1-\frac{W}{L}}

}}}

It can be shown that \xi is not a constant but depends on the parameter !SlopeLong S and the Leakage L. Given that Slope S depends geometrically on Dist D (be more specific), it is not surprising that include Slope improves the calculation of xi. For images with leakage >0 the ratio W/L is deteriorated (Here is room for improvement as the calculation does not take the image orientation into account!). Therefore, xi is calculated as

{{{#!Latex

\xi := \xi_0 + \xi_1\cdot S + \xi_2\cdot\left( 1-\frac{1}{1+\xi_3\cdot L}\right)

}}}

with the following parameters \xi_0=1.299°, \xi_1=0.0632° (°/ns), \xi_2=1.67972° and \xi_3=4.86232.

Further dependencies have been investigated but no significant dependency has been found. In particular, no dependence on the zenith angle could be found.