Issue 
A&A
Volume 618, October 2018



Article Number  A99  
Number of page(s)  11  
Section  Numerical methods and codes  
DOI  https://doi.org/10.1051/00046361/201832927  
Published online  17 October 2018 
Observations of meteors in the Earth’s atmosphere: Reducing data from dedicated doublestation wideangle cameras
^{1} Department of Geodesy and Geoinformation Science, Technische Universität Berlin, Strasse des 17. Juni 135, Berlin, Germany
email: anastasios.margonis@tuberlin.de
^{2} Armagh Observatory, College Hill, Armagh, BT61 9DG UK
email: aac@arm.ac.uk
^{3} Germany Aerospace Center, Institute of Planetary Research, Rutherfordstr. 2, 12489 Berlin, Germany
email: juergen.oberst@dlr.de
Received:
28
February
2018
Accepted:
9
July
2018
Meteoroids entering the Earth’s atmosphere can be observed as meteors, thereby providing useful information on their formation and hence on their parent bodies. We developed a data reduction software package for double station meteor data from the SPOSH camera, which includes event detection, image geometric and radiometric calibration, radiant and speed estimates, trajectory and orbit determination, and meteor light curve recovery. The software package is designed to fully utilise the high photometric quality of SPOSH images. This will facilitate the detection of meteor streams and studies of their trajectories. We have run simulations to assess the performance of the software by estimating the radiants, speeds, and magnitudes of synthetic meteors and comparing them with the a priori values. The estimated uncertainties in radiant location had a zero mean with a median deviation between 0.03^{∘} and 0.11^{∘} for the right ascension and 0.02^{∘} and 0.07^{∘} for the declination. The estimated uncertainties for the speeds had a median deviation between 0.40 and 0.45 km s^{−1}. The brightness of synthetic meteors was estimated to within +0.01 m. We have applied the software package to 177 real meteors acquired by the SPOSH camera. The median propagated uncertainties in geocentric right ascension and declination were found to be of 0.64^{∘} and 0.29^{∘}, while the median propagated error in geocentric speed was 1.21 km s^{−1}.
Key words: meteorites, meteors, meteoroids / techniques: image processing / techniques: photometric / astrometry
© ESO 2018
1. Introduction
Observations of meteors in the Earth’s atmosphere shed light on the properties of the population of meteoroids intercepting the orbit of our planet. The study of the temporal and spatial distribution of meteors requires sensitive optical systems that are able to monitor the night sky. Double station observations (i.e. observations of two cameras from different positions) are required to determine the trajectories and orbit parameters of the meteors.
While algorithms for meteor data reduction are well established in the literature (Ceplecha 1987; TrigoRodríguez et al. 2004; Weryk et al. 2008; Jenniskens et al. 2011), every camera may require an analysis system to account for the specific capabilities of the camera.
For observations of meteors in recent years, our team has used the Smart Panoramic Optical Sensor Head (SPOSH) camera (Oberst et al. 2011). The instrument features a highly sensitive CCD chip that delivers images of high photometric quality. With the wideangle lens, the camera easily captures several hundreds of stars in one image, which requires sophisticated geometric calibration procedures.
To process the data from SPOSH, we developed a comprehensive software package, which allows us to carry out camera calibration, meteor detection, meteor trajectory determination, meteor photometric modelling, and orbit reconstruction. In this paper, we describe this software and demonstrate and assess its performance on synthetic and actual meteor data.
2. The SPOSH camera
The SPOSH camera was designed to image faint transient noctilucent phenomena, such as aurorae, electric discharges, meteors, or impact flashes on dark planetary hemispheres from an orbiting platform (Oberst et al. 2011). The camera is equipped with a highly sensitive backilluminated 1024 × 1024 CCD chip and has a custommade optical system of high lightgathering power with a wide field of view (FOV) of 120 × 120^{∘}. The SPOSH camera system is accompanied by a sophisticated digital processing unit (DPU) designed for realtime image processing and communication with a spacecraft. Owing to the allsky coverage and excellent radiometric and geometric properties of the camera, a large number of meteors can be obtained for reliable event statistics.
For outdoor tests and meteor monitoring, the camera is typically mounted on a tripod pointed vertically up at the sky taking one image every 2 s. For the determination of the meteor velocity, a mechanically rotating shutter with a known frequency is mounted in front of the camera lens. The shutter consists of two blades and has a rotating frequency of 250 RPM resulting in an exposure time of 0.06 s for every shutter opening. Doublestation observations have been carried out routinely providing a large dataset of meteor images.
3. Data reduction
The reduction of the meteor data is performed by different software modules. The calibration software SPOSH Calib is a standalone software for the geometric calibration of SPOSH images (Elgner et al. 2006). The trajectory determination module is based on the MOTS software (Koschny & Diaz del Rio 2002), which was initially modified to process SPOSH data (Maue et al. 2006). All modules were developed anew within the scope of this study. The interaction between the different modules can be seen in Fig. 1.
Fig. 1.
Flow chart showing the different modules of the software package. The camera calibration software is used as a standalone program and in the flowchart is depicted as a rectangle with a white background. 
3.1. Meteor detection
Unlike video cameras, where a meteor only spends a fraction of their trajectory in each frame, exposures longer than one second often capture the whole meteor (e.g. a Perseid) in one image. The meteor detection algorithm that we used is based on the Hough transform technique for extracting linear features within images. This method has been used to detect meteors in photographic image data by previous authors (Trayner et al. 1996; Gural 1997).
The algorithm that we developed first generates 8 bit difference images between three consecutive frames, removing the background and highlighting only short temporal variations. Possible nonrelevant information depicted in the margins of the images (e.g. surrounding mountains and manmade structures), typical within large FOVs, are removed by applying a circular mask. Background noise and stellar scintillation are filtered out by first applying an empirical threshold and then a median filter to the image, thus reducing the overall computation time of the algorithm. Each line, represented by a combination of ρ and θ, passing through each of the remaining pixels contributes to the parameter space H(ρ, θ), known also as voting space, by adding the value A_{xy} = 1,(1)
where ρ is the distance from the origin to the closest point on the straight line, and θ is the angle between the x axis and the line connecting the origin with that closest point. The square brackets [ ] indicate rounding to the nearest integer and the normal representation of a line is(2)
The event detection algorithm is triggered each time a certain threshold value is exceeded. This value is compared against the maximum value found in the parameter space of each image and represents the number of pixels lying on a line in the image space.
Several criteria are used to mitigate the effect of false detections. Slow moving objects, such as airplanes and satellites, appear with a characteristic negativepositivenegative pattern in the difference images (Fig. 2). This pattern is compared against a predefined signal by the user, simulating the path of an airplane projected on the image plane. In this way, events appearing in three consecutive images moving with low apparent angular velocities of 0.6^{∘} > ν _{ang} > 2.2^{∘} are rejected as slowmoving objects. This condition also affects meteors appearing close to their radiant position and/or close to the horizon. For every event, its time of occurrence, the central position of the line, and its direction within the image are saved together with the object name (meteor, slowmoving object, or star) in a text file.
Fig. 2.
Difference image showing a meteor trail and an airplane with its characteristic negativepositivenegative pattern in the lower part of the image. 
A quality parameter q_{md} was introduced to determine the threshold value for the Hough transform. The value of the threshold should ideally detect all meteors in the image data when applied to a meteor detection algorithm. At the same time slowmoving objects and random noise patterns resembling lines should be filtered out. To select a suitable value for this parameter, we balance the number of false positives against the number of meteors the algorithm failed to detect (false negatives) and the processing time it takes for the algorithm to scan the images. The value is computed, after applying various weights to the observed quantities, as follows:(3)
where p_{1} is the percentage of the detected meteors, p_{2} the percentage of false detections, p_{3} the processing time in minutes, and w_{1}, w_{2}, and w_{3} the respective weights. The quality parameter is scaled to values between 0 and 1 (Eq. (3)). The maximum value (q_{md} = 1) for a threshold is reached when all meteors are detected (p_{1} = 100), with no false detections (p_{2} = 0), within a userdefined processing time.
We tested the performance of our algorithm using various threshold values and applying these to > 14 000 images corresponding to eight hours of data from two observing sites. The results were compared with meteors identified after visual inspection of the images. The highest value of the quality parameter for this dataset was found for a threshold of 23 with w_{1}:0.6, w_{2}:0.3, and w_{3}:0.1. The threshold value corresponds to the highest number of pixels lying on a line in a given image. Applying these parameters, 70% of the visually identified meteors were successfully detected by the algorithm (true positives), while 15% of the detected events were false detections (false positives). Figure 3 shows the calculated quality parameter for our dataset. High values are computed from data with a relative high signaltonoise ratio in terms of detected meteors and false detections. This performance of the algorithm can be achieved under favourable weather conditions.
Fig. 3.
Quality parameter values computed from 8 h of image data with respect to different threshold values. The faint lines show the values of the quality parameter q_{m} for each of the 8 datasets while the red line shows the mean value of the quality parameter for the different threshold values. 
3.2. Astrometry
3.2.1. Camera calibration
The geometric calibration of the camera is performed by the SposhCalib software in a semiautomatic process using standard stars in the SPOSH images (Elgner et al. 2006). Stars are ideal calibration targets owing to their high abundance in the images and the precise knowledge of their position at a given time. The SPOSH images may feature up to several thousand stars, which are on average equally distributed over the whole image except image corners. By comparing the actual stars in the image with their expected positions based on a priori information about pointing and interior camera parameters, these parameters can be updated in a leastsquares fashion. This provides an accurate knowledge of the interior, i.e. focal length and geometric distortion, and the exterior orientation (pointing) parameters. The coordinates of the stars are taken from the Tycho2 and HIPPARCOS star catalogues (ESA 1997a).
The transformation equations between the image coordinate system (x,y) and the camera coordinate system (X_{cam}, Y_{cam}, Z_{cam}) are described applying an equidistant camera model (Ray 1994),(4)
A high number of standard stars is achieved by performing initially a precalibration with the help of at least six reference stars selected by the user. The precalibration step provides approximate values for the unknown parameters. After this step, a global calibration is performed using all point sources identified as standard stars in the image.
The SPOSH images show significant radial and nonsymmetrical distortion, mathematically expressed as(5) (6) (7) (8)
The equations above describe the radial (Eq. (5)) and nonsymmetrical distortions (Eq. (6)) and the deviations of the image coordinate system from an orthogonal, uniformly scaled coordinate system (Eqs. (7), (8)). The outer and inner orientation of the camera and the distortion parameters introduced by the lens are determined by fitting a 6thorder polynomial function. These distortion terms are added directly to the pixel coordinates of the stars.
The average residual error for the star positions after the calibration is usually less than 0.25 pixel or 1.68^{′} and usually consistent over the whole image. The displacement Δx_{ij} and Δy_{ij} in image coordinates due to radial distortion are stored in two separate TIFF files. These files serve as lookup tables in the subsequent steps providing the undistorted position of each pixel.
3.2.2. Meteor path on the image plane
The projection of a trajectory of a meteor on an image plane can be seen as a meteor trail. By extending the trajectory before and after the luminous path, a line can be defined on the image plane representing the projection of that extended path. Once the line is defined in both images, its radiant can be determined (Sect. 3.3.1).
In order to speed up the process of defining the meteor line, a threshold is applied to each raw meteor image. The threshold is defined at 2σ of the noise level. The use of a relative low threshold ensures that fainter pixels belonging to the meteor trail are considered in the computation of the line. The line parameters defined in the meteor detection procedure (Sect. 3.1) are used to remove unwanted features in each meteor image, considering the proximity of each pixel to the detected line, its distance to the pixel with the maximum votes in the Hough transform, and its intensity value (Fig. 4).
The positions of the remaining pixels are corrected for radial distortion by retrieving pixeloffset values from the lookup tables generated in the calibration step (Sect. 3.2.1). Owing to the equidistant projection model used by the lens system of the camera, perspective distortions in the image are evident that deflect the path of objects moving along a great circle from a straight line to a curved line. To efficiently detect linear features in the image, pixel coordinates are converted from an equidistant to a gnomonic projection, where straight lines in space preserve their straightness when projected on the image plane.
Fig. 4.
Plots showing the consecutive processing steps of a meteor image for the removal of structures not belonging to the meteor trail. From top left to lower right: detected meteor line and pixel with maximum votes (red) in thresholded image, median filter, computed coefficients for each pixel. High intensities represent meteor pixel and selected meteor pixel. 
The line along which the meteor is moving is computed by applying a customised Hough transform. As an input, we use the corrected image coordinates (in subpixel accuracy) belonging to the meteor trail. Lines running diagonal to the meteor trail results into a higher value in voting space than those parallel to the meteor motion, since more pixels lie along the diagonal line (Fig. 5). We handle this effect as follows: First we apply a Hough transform to determine the top 20 lines intersecting the highest number of pixels. Then we perform a weighted Hough transform considering the intensity values. Unlike the standard Hough transform method, which searches for the line with the maximum votes V in parameter space, we defined a ratio coefficient calculated as the sum of intensity values I with respect to the number of pixels that are pixels apart from each line parameter combination. A distance of pixel is needed to identify which pixels lie on the line since the line does not cross the pixel centre (defined at 0.5 pix),(9)
The bestfitting line is defined as the line with the highest ratio. Since the point spread function (PSF) of an imaging system spreads the light of point sources to neighbouring pixels, the light emitted by a meteor also spreads to pixels located perpendicular to its motion. In order to account for the signal within these pixels, we define a buffer zone of 10 pixels perpendicular to the best line computed.
Occasionally, residual features may be located along the buffer zone. As a result, these remaining pixels affect the determination of the meteor line. In order to remove these features, the consecutive pixeltopixel distances are determined revealing gaps between features. Distances higher than a threshold indicate different pixel entities, where entity is a feature consisting of at least two neighbouring pixels; for example, the meteor trail or meteor segment is such an entity. Assuming that the meteor entity has the maximum number of pixels, we remove all secondary features from the line zone. Finally, the meteor line is determined using weighted least squares. The line parameters (slope plus intercept) and the middle point of the meteor trail, defined as the median of the chosen pixel coordinates, are saved in a text file.
Fig. 5.
Simplified meteor example represented by three intensity levels with the dark grey area corresponding to low dn values. The line intersecting the meteor in the left example has the highest value in voting space, while the right line produces the highest ratio and it is the desired outcome. A buffer zone with a width of 20 pixel parallel to the determined line is depicted by the two parallel thin lines. 
3.2.3. Transformation to the spatial trajectory of the meteor
From the estimated parameters of the meteor line, the underlying image points are generated in subpixel accuracy and transformed from the gnomonic back to an equatorial projection. The pixel coordinates x_{c}, y_{c}, are normalised using the parameters of the interior orientation of the camera, i.e. (Sect. 3.2.1)(10)
where x_{p}, y_{p} are the intersection of the optical axis with the image plane (principal point), p_{x}, p_{y} the pixel size, and f the focal length of the camera system. The points are first projected to the camera coordinate system (∣x ∣ = 1) using the equidistant projection equations. The vectors are then transformed to the local (horizontal) coordinate system,(11)
where R_{ωϕκ} is the 3d rotation matrix that relates the camera to the local coordinate system. Finally, the pointing vectors are transformed to the common Earthcentred, Earthfixed (ECEF) coordinate system. The zaxis becomes parallel to the north pole by rotating the local system by an angle 90 − ϕ_{geo} around the xaxis, with ϕ_{geo} the geocentric latitude of the camera location. The xaxis aligns with the direction of the prime meridian after rotating the system around the zaxis by an angle λ_{geo}, where λ_{geo} the geocentric longitude of the camera location,(12)
3.3. Trajectory determination
3.3.1. Meteor geometry
The trajectory of the meteor is determined using the 3D unit vectors of the defined points on the meteor line. These vectors are generated for each camera from the known camera orientation. The vectors point to the meteor trail and are defined in the geocentric coordinate system. Since the meteor line is initially defined in the images using a gnomonic projection, the intersection points of the direction vectors with a unit sphere lie on a great circle. A plane is fitted through all the unit vectors from each station by solving the standard plane equation using leastsquares(13)
where x_{0} = 0 is the origin of the geocentric coordinate system, x is the direction vector, ⟨n_{x}, n_{y}, n_{z}⟩ are the vector components of the normal vector n, and d is the distance from the plane to the origin and in this equation is equal to zero. The apparent radiant RA_{app}, Dec_{app} of the meteor is calculated as the crossproduct of the two normal vectors n_{1} × n_{2}, determined in (13) with the subscripts indicating the two camera stations. The mean altitude of the meteor is computed by intersecting the direction vector of the central point of the meteor from the shuttered meteor station with the plane generated from the direction vectors to the meteor trail from the second station.
To determine the speed and duration of a meteor, each shuttered meteor image is compared with a database of synthetic meteors (see Sect. 5). These meteors have a fixed geometry and orientation with respect to the camera, i.e. the meteor is moving parallel to the xaxis of the camera system and at 100 km above the camera. The projection of the meteor position at time interval t = dt/2 coincides with the principal point of the camera. The database is created by varying two parameters: the speed and duration. The step size of the database is 0.1 km s^{−1} for the velocity and 0.02 s for the duration of the meteor. From the known geometric relation between the image and meteor plane, the meteor image is transformed so that the meteor plane becomes parallel to the image plane and the distance between principal point and plane is adjusted to 100 km (Fig. 6). This normalised image is then compared with synthetic meteor images in the database accounting for (n_{s/0.1)} × (n_{d/0.02}) different combinations for speed and duration, where n_{s} and n_{d} are the resolution of our partitioning in speed and duration, respectively. For each combination, the meteor trail is timeshifted by 0.06 s to account for various beginning points. The Pearson correlation coefficient is calculated between a synthetic meteor in the database and the normalised image. For the best match we follow a topdown searching approach: first a coarse search is made and then gradually the step size is decreased around the parameters showing a higher correlation. The speed (V_{obs}) and duration of the meteor are derived from the synthetic image with the highest correlation.
Fig. 6.
Left panel: reconstructed meteor plane using the determined meteor orientation and position in camera coordinate system. Right panel: normalized meteor plane being parallel to the image plane and at 100 km distance. The crosses in red color highlight the meteor line on both planes. 
A meteoroid experiences a deceleration when it reaches the denser layers of the Earth’s atmosphere. This effect, socalled atmospheric deceleration, depends on the initial speed of the meteoroid and is more prominent for slower meteoroids. In our studies, we are focussing on the fastmoving Perseid meteoroids and therefore, deceleration is ignored here. The Earth’s rotational velocity contributes an extra 0.004^{∘}/s to the calculated right ascension angle of the radiant and is also neglected in this study.
The speed of a meteoroid slightly increases as soon as it experiences the Earth’s gravitational attraction, a phenomenon known as zenithal attraction. This effect is computed by performing two integrations following Jenniskens et al. (2011): one integration backwards in time including the gravitational effects of the EarthMoon system until the meteoroid reaches the Earth’s sphere of influence and a second integration forwards accounting only for the masses of the Sun and the planets. As input for both integrations the state vector of the meteoroid is used. The new state vector yields the position and velocity of a meteoroid at the time it was recorded in the absence of the EarthMoon system. The velocity vector now points to the geocentric radiant (RA_{geo}, Dec_{geo}).
3.4. Heliocentric orbit
The orbital path of a meteoroid around the Sun, requiring knowledge of Earth and Sun positions, is computed using standard solar system ephemerides (DE421). We use the SPICE software library (Acton et al. 2011) to access ephemeris data and retrieve the following geometric transformations. First the state vectors are transformed from an Earthcentred to a Suncentred ecliptic coordinate system in J2000, i.e.(14)
The heliocentric position and velocity vector are simply computed as the vector sums(15)
where r_{met} and r_{earth} are the state vectors of the meteoroid and the Earth with respect to the Sun in the heliocentric ecliptic coordinate system. Finally, the osculating elements of the orbit are determined using the heliocentric state vector from SPICE routines.
3.5. Photometric reduction
3.5.1. Meteor photometry
Photometric information on meteors is extracted by deconvolving the emitted light of a meteor from the registered signal in equal time intervals (Christou et al. 2015). We remove the effects of radial distortions in the raw image and resample it using inverse distance weighting interpolation. The displacement values Δx_{ij} and Δy_{ij} for each pixel are determined from the geometric camera calibration (Sect. 3.2.1). To speed up the interpolation process, we limit the interpolation to a rectangular area around the meteor trail, for which the position is defined (see Fig. 4). The change of the angular velocity of a meteor owing to perspective distortion is taken into account by projecting the previously determined 3D meteor path (Sect. 3.2.3) to the image.
The number of time intervals n_{t} for which the brightness of the meteor is estimated is computed as the ratio of the length of the rectangular area to the spatial sampling resolution defined by the user. A constant spatial sampling size ensures a stable numerical solution for meteors with low angular velocities, but at the same a highresolution photometric profile for meteor with high angular velocities. From the estimated meteor velocity, the time the meteor needs to cross the rectangular area is calculated following an iterative process. The photometric model can now be applied to the meteor line in the interpolated image.
3.5.2. Photometric calibration
For photometric calibration we use stars depicted in the image. Their positions in the image (pixel coordinates) are computed using the DAOPHOT routines (Stetson 1987) and transformed to the equatorial coordinate system at J2000. The stars are identified by querying the VIZIER database (ESA 1997b) and matching them to the brightest stars (m < 8) found within a radius of 30 arcmin (∼4 pixel) from their position. The flux of each star is measured by defining three circular areas around the light source: an inner circular area measuring the light coming from the star and an outer ring determined by two circular areas defining the sky background. We set the star aperture to a radius of 3 × FWHM, which encloses nearly 100% of the stellar flux (Merline & Howell 1995). The instrumental magnitude m_{inst} is then defined as(16)
where A is an arbitrary constant, C_{i} is the DN value in the i^{th} pixel, C_{sky} is the mean sky background value, n is the number of pixels in each aperture, and t is the integration time of the frames.
To transform the computed instrumental magnitudes to a standard photometric system, we first convert the HIPPARCHOS H_{p} magnitudes from the Vizier database into Johnson V magnitudes using the following expression (Harmanec 1998):(17)
where B − V is the colour index of each star from the VIZIER database and α_{i} the transformation coefficients. The light emitted by each star is partially absorbed by the Earth’s atmosphere. Therefore, the amount of the absorption for each star is proportional to the amount of atmosphere the light has to traverse to reach an observer on the Earth’s surface. This means that light of stars appearing close to the horizon experiences a greater absorption than stars close to the zenith. The amount of atmosphere, called airmass, is calculated as(18)
where z is the zenith angle (Binzel 2006). The equation is taking into account the curvature of the Earth. Moreover, the attenuation of light is computed as a function of the wavelength due to Rayleigh scattering and therefore, the amount of attenuation for each star depends on its colour. To account for the colour difference of our star field, we apply a colour correction for each star using the colour indices from the star catalogue. The instrumental magnitudes are corrected for atmospheric effects and converted to absolute magnitudes,(19)
where X is the airmass, T_{c} the transformation coefficient, CI the colour index, k is the extinction coefficient given in magnitudes per unit airmass, and Z_{p} is a scaling factor. We compute three sets of correction parameters for U, V, and I colour corrections using least squares.
We tested our photometric calibration module with a typical SPOSH image. We detected 393 stars in the image, where the faintest is of +6.3 mag. A high correlation between calibrated and standard stellar magnitudes using the V − I colour index was found, matching the spectral response of the system. Stars with an airmass greater than 3 were excluded from the procedure. Figure 7 shows the calibrated magnitudes of the stars with respect to their catalogue magnitudes from a single frame. The standard deviation between catalogue and measured magnitudes was of 0.22 mag.
Fig. 7. Calibrated star magnitudes vs. standard star V magnitudes. Dashed line shows the ideal onetoone relationship between the two quantities. 
4. Error propagation
The errors of the unknown parameters are calculated by applying error propagation. In the sections to follow we refer to these as the propagated uncertainty (or propagated errors) to distinguish this uncertainty from the statistics of differences between estimated and a priori known parameters (estimated uncertainty) as well as the uncertainty in estimating a common property of the meteors, for example the radiant and speed of a shower, by taking the average over a number of meteors (observed uncertainty). We encounter the estimated uncertainty principally in tests with our synthetic data (Sect. 5.1). The unknown parameters are the apparent and geocentric radiant positions, observed, geocentric, and heliocentric speed of the meteoroid, and orbital elements of its orbit around the Sun. The observed quantities are the parameters of the meteor line ρ and θ. The general law of error propagation is of the form(20)
where C_{xx} is the stochastic model of the measurements and y are the parameters to be estimated. The parameter C_{yy} is the variancecovariance matrix of the unknown parameters. The uncertainties of the direction and location of the meteor line on the image affect the uncertainties of the parameters and need to be carefully estimated. We used our synthetic meteor dataset (see Sect. 5) to estimate the line uncertainties. The distribution of the estimated uncertainties for ρ and θ are well approximated by a normal distribution with a standard deviation of 0.07^{∘} for θ and 1.35 pixel for the distance of the projection to the line. These values depend highly on the length of the meteor trail, PSF, and resolution of the CCD.
5. Software validation
5.1. Synthetic meteor data
We verified our software modules with the help of synthetic meteor data. A meteor trail is generated by providing a number of parameters which i) define the dynamic and photometric properties of a meteor, ii) define the geometric relation between the observers and the meteor, and iii) projects the luminous path of the meteor to the image plane of the given camera system. Table 1 summarises the initial conditions used to generate synthetic meteor trails.
Once the meteor path is generated in space, the corresponding meteor trail is projected on the image plane. The trail is created by convolving a 2D Gaussian curve imitating the motion of a pointlike light source on a given camera system. The method is based on the photometric model in Christou et al. (2015) implemented in reverse. The peak intensity value of each meteor is kept constant while the standard deviation of the Gaussian PSF is set equal to one pixel. The brightness of each synthetic meteor is normally distributed along the meteor trail. The peak brightness also varies between each meteor and resembles different shape curves (Beech & Hargrove 2004; Borovička et al. 2007). The position of the peak along the meteor trail in our sample follows a normal distribution with its centre at n_{t}/2 and a standard deviation of n_{t}/10, where n_{t} is the number of time intervals. The meteor trails in one of the stations were chopped periodically to simulate the effect of the rotating shutter. The starting point of a meteor at time t_{0} is placed randomly within a shutter break and ranges between 0 s and 0.06 s. As an example, a meteor with t_{0} = 0 will receive light directly, while a meteor with t_{0} = 0.06 occurs exactly at the time when the shutter is located in front of the lens. For the parameters of the inner orientation of the camera, typical values for the SPOSH camera were used. The pointing of both cameras was chosen to slightly deviate from an optical axis parallel to the zenith. The distance between the two stations was set at 55.6 km. All synthetic meteors in our simulations occur at the same time, i.e. time information is not relevant. The position of the radiant is given in equatorial coordinates system as RA and Dec.
An image with random noise was generated using the noise properties of the SPOSH images that are the same size as the meteor image. Additionally, 20 2DGaussian PSF simulating the stellar sources in the image were distributed randomly in the FOV of the camera. The image with the synthetic stars was added to the noise image. The position and brightness levels of the stars were kept fixed for all synthetic images. The light of the synthetic stars in the noise image represent the remaining light due to scintillations, visible in the SPOSH difference image data. Finally, the noise image was added to the meteor image, creating the input for our algorithm. The software uses the images of each synthetic meteor as input and calculates its radiant position, speed, brightness, and heliocentric orbit.
We generated a dataset of synthetic meteor trails considering different geometric configurations. We present results for two types of synthetic data. For the first, we created 208 synthetic meteors with random positions and directions with respect to the location of the cameras, and then used our program to estimated their radiants, speeds, and magnitudes. Fortyseven of the synthetic meteors had convergence angles Q ≤ 10^{∘} yielding large errors in the radiant determination. These meteors were therefore excluded from the procedure. For the remaining 161 synthetic meteor pairs with Q > 10^{∘} we determined the radiant position (RA and Dec) and the speed (V) and computed their estimated uncertainties as the difference between the a priori value and that calculated from the code. We describe the statistical dispersion of the probability distributions by calculating the median absolute deviation (MAD), which is statistically a more robust measure for asymmetric distributions than the standard deviation. The distributions of the estimated uncertainties for RA and Dec were centred at zero with a median deviation of 0.11^{∘} and 0.07^{∘}, respectively (Fig. 8). The statistical properties of the estimated uncertainties are shown graphically in Fig. 9. The propagated uncertainties had a median value of 0.33^{∘} for RA and 0.16^{∘} for Dec (Fig. 10). The propagated and estimated errors are in good agreement, which implies a realistic stochastic model (Table 2). The distribution of the estimated errors for the speed appears to be offset from zero with a median of 0.24 km s^{−1} and a median deviation of 0.51 km s^{−1}. For calculating the meteor magnitudes, we used a subset of the data consisting of 185 unshuttered meteors of which the individual residual is ≤0.5 m. The estimated uncertainties had a median value of +0.01 m and a median deviation of 0.03 m.
Fig. 8.
Left panel: radiant dispersion for synthetic meteors originating from random directions. Right panel: radiant dispersion for meteors with the radiant point located in the local zenith. The filled circles (•) represent meteors appearing 50^{∘} above the horizon while open circles (∘) show meteors with elevation angles lower than 50^{∘}. 
Fig. 9. Distribution of the estimated uncertainties in RA, Dec and speed for the synthetic meteors. The length of the boxes indicates the dispersion of the data. Each box encloses 50% of the data. The extending vertical lines from the boxes indicate the range of 80% of the data with the lower and upper horizontal bars marking the 10% and 90% levels. Data outside the 80% range are shown as open circles (∘). The horizontal line inside each box indicates the median value. The units are degrees for RA and Dec, and km s^{−1} for speed. 
Fig. 10. Distribution of propagated uncertainties for the synthetic meteors. For a description of the plot see the caption of Fig. 9. 
A second set of 208 synthetic meteors was then created with the same radiant point for all meteors placed at the local zenith to simulate a meteor shower observed by the two cameras. One hundred seventy of the meteors had a convergence angle Q > 10^{∘}. Nine radiants with large estimated uncertainties were excluded from the procedure. As for the first set of synthetic data, the estimated errors for RA and Dec also have a zero median but slightly lower median deviation of 0.05^{∘} and 0.03^{∘}, respectively (Fig. 9). These reduce to 0.03^{∘} and 0.03^{∘} when considering only 126 meteors occurring > 40^{∘} from the local horizon. The dispersion of the propagated uncertainties for RA and Dec are similar to those computed for the first synthetic dataset. The median propagated uncertainty was 0.37^{∘} for RA and 0.27^{∘} for Dec (Fig. 10). Figure 11 shows the relation between the position of a meteor in the image and the estimated errors in radiant. Since the pointing of the camera and the radiant of the simulated shower was set to the local zenith, the angular separation between the radiant and position of a meteor in the image corresponds to the elevation angle of the meteor. The median for the estimated errors in the speed for these 161 meteors was 0.12 km s^{−1}. The median of the estimated errors for 192 unshuttered synthetic meteors, for which the residuals in magnitude Δmag were < 0.5 m, was +0.01 m with a median deviation of 0.02 m.
Fig. 11. Plot showing the relation between the angular separation and the estimated uncertainties in right ascension (filled circles) and declination (open circles). 
Different types of uncertainties for the synthetic and real meteor data.
5.2. Real meteor data
We applied our software to three hours of doublestation SPOSH image data acquired during an observing campaign held in Greece (i.e. at ∼37^{∘}N latitude ) on 12 August 2015 from 23 to 02 UT. The cameras were pointing to the zenith with their xaxis orientated to the north. The baseline between the two sites was 51.5 km. We reduced 177 meteor image pairs and determined their trajectories, velocities and heliocentric orbits.
We focus on a 20 × 20^{∘} area centred at RA = 46^{∘}, Dec = 58^{∘}, close to the nominal radiant position of the Perseids (Fig. 12). To distinguish between Perseid and nonPerseid meteors, we performed a classification based on radiant position and speed as follows: 132 meteors were found to radiate from within this area. We determine the radiant of the Perseid shower as the median value of these radiants: RA = 45.96^{∘} and Dec = 57.77^{∘}. We assume that most of the meteors are Perseids and we find the 1σ uncertainty in RA and Dec to be 3.29^{∘} and 2.27^{∘}, respectively. The geocentric speeds V_{G} have a median value at 58.97 km s^{−1} and a median propagated error of 1.21 km s^{−1} (Fig. 13). We classify all meteors with speeds closer to this median speed than four times the median propagated error, as Perseids. In this way, we identified 71 meteors belonging to the Perseids meteor shower (Fig. 12). Their median speed V_{G} was found to be 59.58 km s^{−1} with a median deviation of the observed uncertainties of 0.48 km s^{−1}. Statistical properties are given in Table 2. As the aim of our example is to demonstrate the successful usage of our software using real data, we neglected the effect of radiant drift.
Fig. 12. Geocentric radiants of 132 meteors originating close to the Perseid radiant. The ellipse defines the area occupied by Perseids. Meteors classified as Perseids according to their speed are shown as filled circles (•). All other meteors not meeting the requirements to be classified as Perseids are shown as open circles (∘). 
Fig. 13. Speed distribution of 132 meteors originating from an area close to the Perseid radiant. The dashed line shows computed median speed of 59.58 km s^{−1}. The bars in dark grey show the distribution of speeds for 71 meteors classified as Perseids found within the same area. 
We calculated the magnitudes of the 71 meteors identified as Perseids from the unshuttered images. We defined this magnitude to be the brightest value obtained for the light curve of each meteor. The magnitude distribution index r for the shower was found to be 2.10 ± 0.10 (Fig. 14). The mean value is slightly above the upper limit of the range 1.87 < r < 2.01 given by Brown & Rendtel (1996) during the years 1988–1994 but within the range 1.86 < r < 2.12 found by Jenniskens et al. (2016) for the years 1991–1994. These observations cover both the pre and post activity of the meteor shower that coincided with the perihelion passage of the parent comet in 1992. The index was computed from 61 meteors brighter than +0 m. The light curve of a bright, doubleflaring Perseid was computed following the method described in Sect. 3.5.1 (Fig. 15). The time step size dt was 0.006 and 0.01 s for the unshuttered and shuttered camera respectively. The gaps between the points represent the time intervals with no information owing to the rotating shutter.
Fig. 14.
Left panel: magnitude distribution of all meteors located inside the ellipse in Fig. 12. Right panel: cumulative distribution of the magnitudes. The slope of the straight line defines the mass distribution index r for the shower during the observing time. Only meteors brighter than +0 m were included in the fit. The dashed line represents the cutoff value. 
Fig. 15.
Brightness profile as a function of time for a bright Perseid meteor. The line shows the absolute magnitude calculated from the continuous meteor trail while the red dots represent brightness calculation from the shuttered meteor trail. We note the two flares at ∼0.9 s and ∼1.0 s. 
Heliocentric orbits were computed for all 177 meteors in our sample. The median orbital elements of the 71 Perseid meteors were compared with the orbits found in five studies (Kresák & Porubčan 1970; Jopek et al. 2003; SonotaCo 2009; Jenniskens et al. 2016) and the orbit of comet 109P/SwiftTutle, parent comet of the Perseids (Jenniskens 2006) (Table 3). In general, we find a good agreement in the radiant, speed, and orbital elements. One exception is in the semimajor axis, which is known to be very sensitive to variations in the velocity of a meteoroid (Williams 1996; Jenniskens 1998).
Radiant positions, speeds, and orbital elements of Perseid meteors found in 5 studies compared with median values computed in this work and the orbit of the parent comet.
6. Conclusions
We have presented SPOSHRed, a software package for the data reduction of doublestation meteor image data acquired by the SPOSH camera. The software extracts information about trajectories, heliocentric orbits and brightness levels of recorded meteors.
The software was tested for simulated and real meteor data. We simulated different geometric configurations between a meteor shower and two observing sites. We suggest that such simulations can be used to assess the quality of the derived meteor trajectories for different camera network configurations and predicted meteor shower or outbursts. We expect that the results will greatly contribute to the planning of observing campaigns by finding the best location and orientation between the camera stations and predicted radiant position.
The software presented in this paper was developed to reduce data acquired by the SPOSH camera system. In the future we plan to provide a more generic version of the software package that can handle image datasets recorded by different camera systems.
The real meteor data used in this work is part of a large dataset comprising eight years of Perseids observations using SPOSH. The accumulated observing period spans more than 20 days within the activity period of the Perseids. During this period > 15 000 single meteors have been recorded from both stations. The reduction software opens the opportunity to analyse the available unique SPOSH meteor data. A full analysis of the data focussing on the Perseid meteor shower will be presented in a forthcoming paper.
Acknowledgments
The authors would like to thank the reviewer for providing us with useful comments and detailed suggestions for the manuscript.
References
 Acton, C., Bachman, N., Del Diaz Rio, J. 2011, in EPSCDPS Joint Meeting, 2011, 32 [NASA ADS] [Google Scholar]
 Beech, M., & Hargrove, M. 2004, Earth Moon Planet, 95, 389 [NASA ADS] [CrossRef] [Google Scholar]
 Binzel, R. P. 2006, Minor Planet Bull., 33, 91 [NASA ADS] [Google Scholar]
 Borovička, J., Spurný, P., & Koten, P. 2007, A&A, 473, 661 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Brown, P., & Rendtel, J. 1996, Icarus, 124, 414 [NASA ADS] [CrossRef] [Google Scholar]
 Ceplecha, Z. 1987, Bull. astr. Inst. Czechosl., 38, 222 [Google Scholar]
 Christou, A. A., Margonis, A., & Oberst, J. 2015, A&A, 581, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Elgner, S., Oberst, J., Flohrer, J., & Albertz, J. 2006, Fifth International Symposium TurkishGerman Joint Geodetic Days, 29–31 March 2006, Berlin [Google Scholar]
 ESA 1997a, The HIPPARCOS and Tycho Catalogues, Astrometric and Photometric Star Catalogues Derived from the ESA HIPPARCOS Space Astrometry Mission, ESA SP, 1200 [Google Scholar]
 ESA 1997b, VizieR Online Data Catalog: I/239, [Google Scholar]
 Gural, P. S. 1997, WGN, J. Int. Meteor Organ., 25, 136 [NASA ADS] [Google Scholar]
 Harmanec, P. 1998, A&A, 335, 173 [NASA ADS] [Google Scholar]
 Jenniskens, P. 1998, Earth Planet. Space, 50, 555 [NASA ADS] [CrossRef] [Google Scholar]
 Jenniskens, P. 2006, Meteor Showers and their Parent Comets (Cambridge: University Press) [CrossRef] [Google Scholar]
 Jenniskens, P., Gural, P. S., Dynneson, L., et al. 2011, Icarus, 216, 40 [NASA ADS] [CrossRef] [Google Scholar]
 Jenniskens, P., Nénon, Q., Albers, J., et al. 2016, Icarus, 266, 331 [NASA ADS] [CrossRef] [Google Scholar]
 Jopek, T. J., Valsecchi, G. B., & Froeschlé, C. 2003, MNRAS, 344, 665 [NASA ADS] [CrossRef] [Google Scholar]
 Koschny, D., & Diaz del Rio, J. 2002, WGN, J. Int. Meteor Organ., 30, 87 [Google Scholar]
 Kresák, L., & Porubčan, V. 1970, Bull. Astr. Inst. Czechosl., 21, 153 [NASA ADS] [Google Scholar]
 Maue, T., Flohrer, J., & Oberst, J. 2006, The Meteor Orbit Determination (MOD) Workshop, 11–13 September 2006, Roden, The Netherlands, eds. D. Koschny, & J. McAuliffe [Google Scholar]
 Merline, W. J., & Howell, S. B. 1995, Exp. Astron., 6, 163 [NASA ADS] [CrossRef] [Google Scholar]
 Oberst, J., Flohrer, J., & Elgner, S. 2011, 59, 1 [Google Scholar]
 Ray, S. 1994, Applied Photographic Optics (Oxford, UK: Focal Press) [Google Scholar]
 SonotaCo. 2009, WGN, J. Int. Meteor Organ., 37, 55 [Google Scholar]
 Stetson, P. B. 1987, PASP, 99, 191 [NASA ADS] [CrossRef] [Google Scholar]
 Trayner, C., Haynes, B. R., & Bailey, N. J. 1996, in 3. IEEE International Conference on Image Processing, 2, 693 [CrossRef] [Google Scholar]
 TrigoRodríguez, J. M., CastroTirado, A. J., Llorca, J., et al. 2004, Earth Moon Planet., 95, 553 [NASA ADS] [CrossRef] [Google Scholar]
 Weryk, R. J., Brown, P. G., Domokos, A., et al. 2008, Earth Moon Planet., 102, 241 [NASA ADS] [CrossRef] [Google Scholar]
 Williams, I. P. 1996, Earth Moon Planet., 72, 321 [NASA ADS] [CrossRef] [Google Scholar]
All Tables
Radiant positions, speeds, and orbital elements of Perseid meteors found in 5 studies compared with median values computed in this work and the orbit of the parent comet.
All Figures
Fig. 1.
Flow chart showing the different modules of the software package. The camera calibration software is used as a standalone program and in the flowchart is depicted as a rectangle with a white background. 

In the text 
Fig. 2.
Difference image showing a meteor trail and an airplane with its characteristic negativepositivenegative pattern in the lower part of the image. 

In the text 
Fig. 3.
Quality parameter values computed from 8 h of image data with respect to different threshold values. The faint lines show the values of the quality parameter q_{m} for each of the 8 datasets while the red line shows the mean value of the quality parameter for the different threshold values. 

In the text 
Fig. 4.
Plots showing the consecutive processing steps of a meteor image for the removal of structures not belonging to the meteor trail. From top left to lower right: detected meteor line and pixel with maximum votes (red) in thresholded image, median filter, computed coefficients for each pixel. High intensities represent meteor pixel and selected meteor pixel. 

In the text 
Fig. 5.
Simplified meteor example represented by three intensity levels with the dark grey area corresponding to low dn values. The line intersecting the meteor in the left example has the highest value in voting space, while the right line produces the highest ratio and it is the desired outcome. A buffer zone with a width of 20 pixel parallel to the determined line is depicted by the two parallel thin lines. 

In the text 
Fig. 6.
Left panel: reconstructed meteor plane using the determined meteor orientation and position in camera coordinate system. Right panel: normalized meteor plane being parallel to the image plane and at 100 km distance. The crosses in red color highlight the meteor line on both planes. 

In the text 
Fig. 7. Calibrated star magnitudes vs. standard star V magnitudes. Dashed line shows the ideal onetoone relationship between the two quantities. 

In the text 
Fig. 8.
Left panel: radiant dispersion for synthetic meteors originating from random directions. Right panel: radiant dispersion for meteors with the radiant point located in the local zenith. The filled circles (•) represent meteors appearing 50^{∘} above the horizon while open circles (∘) show meteors with elevation angles lower than 50^{∘}. 

In the text 
Fig. 9. Distribution of the estimated uncertainties in RA, Dec and speed for the synthetic meteors. The length of the boxes indicates the dispersion of the data. Each box encloses 50% of the data. The extending vertical lines from the boxes indicate the range of 80% of the data with the lower and upper horizontal bars marking the 10% and 90% levels. Data outside the 80% range are shown as open circles (∘). The horizontal line inside each box indicates the median value. The units are degrees for RA and Dec, and km s^{−1} for speed. 

In the text 
Fig. 10. Distribution of propagated uncertainties for the synthetic meteors. For a description of the plot see the caption of Fig. 9. 

In the text 
Fig. 11. Plot showing the relation between the angular separation and the estimated uncertainties in right ascension (filled circles) and declination (open circles). 

In the text 
Fig. 12. Geocentric radiants of 132 meteors originating close to the Perseid radiant. The ellipse defines the area occupied by Perseids. Meteors classified as Perseids according to their speed are shown as filled circles (•). All other meteors not meeting the requirements to be classified as Perseids are shown as open circles (∘). 

In the text 
Fig. 13. Speed distribution of 132 meteors originating from an area close to the Perseid radiant. The dashed line shows computed median speed of 59.58 km s^{−1}. The bars in dark grey show the distribution of speeds for 71 meteors classified as Perseids found within the same area. 

In the text 
Fig. 14.
Left panel: magnitude distribution of all meteors located inside the ellipse in Fig. 12. Right panel: cumulative distribution of the magnitudes. The slope of the straight line defines the mass distribution index r for the shower during the observing time. Only meteors brighter than +0 m were included in the fit. The dashed line represents the cutoff value. 

In the text 
Fig. 15.
Brightness profile as a function of time for a bright Perseid meteor. The line shows the absolute magnitude calculated from the continuous meteor trail while the red dots represent brightness calculation from the shuttered meteor trail. We note the two flares at ∼0.9 s and ∼1.0 s. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.