Apparatus, methods, and programs for satellite monitoring

Abstract

A process for satellite tracking including obtaining image data representing a field of view of an optical system, the image data representing one or more satellite trails and a celestial background, and being mapped to a coordinate system in fixed spatial relation to the celestial background; executing differential image processing on the obtained image data to extract one or more satellite trails from the celestial background; and executing matched filtering to calculate orbital parameters representing each of the one or more extracted satellite trails.

Description

[0001] TITLE

[0002] Apparatus, Methods, and Programs for Satellite Monitoring

[0003] TECHNICAL FIELD

[0004] The invention is in the field of satellite technology and in particular relates to calculating orbital parameters of Earth-orbiting satellites based on images.

[0005] BACKGROUND

[0006] Space-based infrastructure plays a crucial role in modem society, providing numerous indispensable services, including satellite navigation and positioning, internet access, disaster management, weather forecasting, forestry and agriculture, as well as scientific research and space exploration. These services are not only integral to the modem and future economy but are also vital for tackling climate change and expanding human horizons. The majority of satellites are deployed in Low Earth Orbit (LEO) due to its proximity to Earth, which is of great importance for both telecommunications and Earth observation applications. However, the LEO region is becoming increasingly congested with defunct satellites and space debris, an issue that also extends to higher orbits. Earth’s orbital environment is a finite resource, necessitating careful planning and close international coordination, as promoted by the United Nations’ guidelines for the Long-Term Sustainability of Outer Space Activities. Accordingly, to improve orbital planning of such a finite resource, it is desirable to achieve a comprehensive global monitoring system of the space environment, including the tracking of active satellites and space debris. This is particularly desirable in the light of the rapid increase in space objects since 2020, coinciding with the rise of reusable delivery systems, as highlighted by the European Space Agency report (European Space Agency, “Space environment report,” European Space Agency, Tech. Rep., 2023.)

[0007] An important aspect of improving satellite safety is the implementation of advanced tracking and monitoring systems to detect and predict the movement of space debris and other objects, with a view to reducing collision risk. This monitoring and tracking is typically referred to as Space Situational Awareness (SSA), and is largely achieved through the use of ground-based sensors. A key SSA sensing mechanism is radar technology, which typically involves transmitting high-power modulated radio pulses and then collecting the reflections from space objects. A return signal arrives after a certain delay and incurs a Doppler shift due to the relative motion. Based on the delay and Doppler shift, radar can simultaneously estimate the distance and radial velocity of the space object. Furthermore, recent advances in software-defined radio chipsets have enabled the development of radars with multibeam capabilities, allowing the detection and tracking of multiple objects simultaneously. However, because the eneigy wavefront in both the transmit and receive paths expands according to the inherent inverse power law, the resulting signal loss is proportional to ρ4, where ρis the distance to the target.Another important method for monitoring satellites is the use of optical sensors, typically performed with the aid of a telescope equipped with a precise azimuth / elevation tracking mount. Unlike radar sensors, the intensity of the captured signal in passive optical sensors is proportional to ρ2since it relies of natural illumination by the sun. However, this method is subject to illumination availability and weather conditions when performed from the ground. Also, passive optical sensors do not provide direct range information; instead, the satellite orbit must be estimated based on observing its angular coordinates across time. Other limitations of optical systems relate to sky brightness and tropospheric conditions, which affect the ability to observe during daytime and in poor weather. As such, space-based (in-orbit) optical systems have advantages with respect to their ground-based counterpart because they do not suffer from atmospheric issues.

[0008] Space-based SSA have long been utilized for space situational awareness, such as the Pathfinder satellite in 2010, equipped with a gimbal camera, and the Sapphire satellite in 2013. Other SSA approaches are based on characterizing the unique Doppler profile of an active RF beacon transmitted by approaching satellites (A. Al-Hourani, " In-Orbit Space Situational Awareness Using Doppler Frequency Shift," IEEE Transactions on Aerospace and Electronic Systems, vol. 60, no. 5, pp. 7542-7547, Oct. 2024,).

[0009] Practical optical systems dedicated to SSA typically have a limited field of view (FoV) to improve angular resolution, where a narrower FoV provides higher resolution and vice versa. Despite the advantages of optical monitoring, this inherent limitation results in practical FoVs on the order of less than a few degrees (T. Yanagisawa, H. Kurosaki, H. Oda et al., “Ground-based optical observation system for LEO objects,” Advances in Space Research, vol. 56, no. 3, pp. 548-564, 2015.). Consequently, the surveillance capability is limited, confining typical optical systems to tracking predetected objects and requiring scheduled tracking tasks (P. M. Siew and R. Linares, “Optimal tasking of ground-based sensors for space situational awareness using deep reinforcement learning,” Sensors, vol. 22, no. 20, 2022. [Online]. Available: >

[0010]

[0011] In summary, shortcomings in existing technologies include:

[0012] Dependency on prior knowledge of orbital parameters;

[0013] High-cost: this is particularly true of satellite-based systems;

[0014] Limited coverage: many existing technologies are limited to narrow FoV sensors, resulting in incomplete or slow tracking capabilities;

[0015] Complex Infrastructure: many existing technologies leverage ground-based radar systems requiring large, costly infrastructure.It is desirable to at least partially ameliorate the constraints of the state of the art, which rely upon narrow FoV imaging. It is desirable to develop robust detection methods for satellite trails in sensors regardless of field-of-view, that is, detection methods that are not limited to narrow FoV optical sensors.

[0016] STATEMENTS

[0017] Embodiments may include a method comprising: obtaining image data representing a field of view of an optical system, the image data representing one or more satellite trails and a celestial background, and being mapped to a coordinate system in fixed spatial relation to the celestial background; executing differential image processing on the obtained image data to extract one or more satellite trails from the celestial background; executing matched filtering to calculate orbital parameters representing each of the one or more extracted satellite trails.

[0018] The method may be a computer-implemented method and may be a method for satellite tracking, satellite monitoring, satellite orbit calculation, satellite orbit estimation, or satellite orbit computation. The method may be implemented by an apparatus comprising memory hardware storing processing instructions and processor hardware configured to execute the stored processing instructions, the execution of the stored processing instructions causing the processor hardware to perform the method. The image data may be image data comprising images of visible light sensed by an optical sensor, or may be images comprising images of infrared light or near-infrared light sensed by an optical sensor.

[0019] The image data may be image data comprising a series of images, each image being a representation of visible light captured by an optical sensor during an exposure of the optical sensor to light from the field of view via optics. The image data be image data comprising a sequence of difference images generated by an optical system such as an event-based camera, which may be referred to as an event camera, a neuromorphic camera, a silicon retina, or a dynamic vision sensor, wherein rather than each pixel value representing intensity of light captured during an exposure, each pixel value represents brightness change at the pixel during an imaging period.

[0020] Embodiments provide for satellite detection and satellite trail calculation that is compatible with narrow FoV optical systems, wide FoV optical systems, and hybrid wide- and narrow- FoV optical systems.

[0021] Embodiments provide apparatus and methods for empirical estimation (i.e. calculation based on observation) of satellite orbits using optical sensor hardware.Embodiments provide a cost-effective, ground-based or space-based, wide-area optical system for realtime satellite detection and orbital estimation. Embodiments operate without prior knowledge of orbits, providing accurate tracking irrespective of whether the optical system is ground-based or space-based.

[0022] Optionally, the calculated orbital parameters comprise values for one or more of the six Keplerian elements. Optionally, the calculated orbital parameters for each of the six Keplerian elements.

[0023] Optionally, executing matched filtering comprises estimating a set of orbital parameter values for the one or more satellite trails in the image data and measuring a discrepancy between the obtained image data and a set of emulated satellite trails generated by the estimated set of orbital parameter, and optimising the estimated set of orbital parameter values through iteration to satisfy a predefined optimisation condition.

[0024] Optionally, executing matched filtering comprises estimating a set of orbital parameter values for the one or more satellite trails in the image data and measuring a discrepancy between the obtained image data and an emulated satellite trail generated by the estimated set of orbital parameter values, and iteratively repeating the estimating and measuring to optimise the estimated set of orbital parameter values to satisfy a predefined optimisation condition.

[0025] Optionally, executing matched filtering comprises estimating a set of orbital parameter values for the one or more satellite trails in the image data and measuring a discrepancy between the obtained image data and a set of emulated satellite trails generated by the estimated set of orbital parameter values and passed through a model representing an optical system via which the image data is generated, and optimising the estimated set of orbital parameter values through iteration to satisfy a predefined optimisation condition.

[0026] Optionally, the iteratively estimating and measuring to optimise the estimated set of orbital parameter values is performed by an optimisation algorithm. The optimisation algorithm may be a particle swarm optimisation algorithm.

[0027] Optionally, the emulated satellite trail is passed through a model representing an optical system via which the image data is generated before the discrepancy is measured.

[0028] Optionally, as each estimated set of orbital parameter values is optimised the corresponding orbit is eliminated from the obtained image data for optimising subsequent set(s) of orbital parameter values.Optionally, the proposed values for the Keplerian elements is a matrix of values comprising, for each of the one or more extracted satellite trails, a matrix element per orbital parameter.

[0029] Optionally, for cases in which there are multiple satellite trails in the obtained image data, the executing the matched filtering comprises:

[0030] (i) estimating a set of orbital parameter values for the first satellite trail in the image data and measuring a discrepancy between the obtained image data and an emulated satellite trail generated by the estimated set of orbital parameter values, and iteratively repeating the estimating and measuring to optimise the estimated set of orbital parameter values to satisfy a predefined optimisation condition;

[0031] (ii) simulating a satellite trail of the optimised set of orbital parameter values in the image data;

[0032] (iii) supressing or eliminating the simulated trail from the image data;

[0033] (iv) repeating steps (i), (ii), and (ii) for each remaining satellite trail. And optionally, (v) outputting the optimised sets of orbital parameters and / or the simulated satellite trails.

[0034] Optionally, the obtained image data comprises images captured by event-based camera hardware (which may be referred to as an event camera, a neuromorphic camera, a silicon retina, or a dynamic vision sensor), the images being a series of difference images representing differences between a series of two or more exposures, and the differential image processing comprises computing a background image from the obtained series of difference images by calculating a mean image of the obtained difference images. The event-based camera hardware may be sensitive to visible light (for example, light with wavelength in the range 380 to 700nm), infrared light (for example, light with wavelength in the range 700nm to 1mm), or near-infrared light (for example, light with wavelength in the range 750nm to 1.4microns). Such image data may comprise readings only for pixels sensing a change in brightness or intensity of received light during an imaging period, with the readings quantifying the change.

[0035] Optionally, the obtained image data comprises a sequence of exposures by optical sensor hardware of the optical system, and executing the differential image processing comprises:

[0036] (i) computing a series of difference images comprising a difference image between each pair of adjacent exposures in the sequence of exposures;

[0037] (ii) computing a background image from the series of difference images by calculating a mean image of the series of difference images. The optical sensor hardware may be an optical sensor configured to sense intensity of light in a defined wavelength received by each pixel of the sensor during an exposure. The optical sensor may be sensitive to visible light (for example, light with wavelength in the range 380 to 700nm), infrared light (for example, light with wavelength in the range 700nm to 1mm), or nearinfrared light (for example, light with wavelength in the range 750nm to 1.4microns).Optionally, the series or sequence of difference images are processed prior to the computing the background image by normalizing the series of difference images and / or bounding the pixel values to a predetermined saturation level.

[0038] Optionally, obtaining the image data comprises, at an optical system, taking one or a series of exposures of a field-of-view of the optical system by reading data from optical sensor hardware.

[0039] Optionally, the image data comprises, for each image, pixel values for an arrangement of pixels representing intensity of light received at each pixel over an exposure period.

[0040] Optionally, the image data comprises a sequence of difference images generated by an optical system such as an event-based camera, which may be referred to as an event camera, a neuromorphic camera, a silicon retina, or a dynamic vision sensor, wherein rather than each pixel value representing intensity of light captured during an exposure, each pixel value represents brightness change at the pixel during an imaging period.

[0041] Optionally, the optical system is satellite-mounted, is airborne, or is spaceborne.

[0042] Optionally, the optical system is ground-based.

[0043] Optionally, the optical system is mounted on an equatorial mount.

[0044] Optionally, the optical system comprises an array of cameras, each camera comprising an optical lens and an optical sensor onto which incident light is focussed by the optical lens.

[0045] Optionally, the optical system comprises a single camera comprising an optical lens and an optical sensor onto which light is focussed by the optical lens.

[0046] Optionally, the optical lens of the single camera is a fish-eye lens.

[0047] The method may further comprise using the calculated orbits for collision avoidance. In such cases, the method may be a satellite control method, a satellite collision avoidance method, or a satellite manoeuvre control method.

[0048] Optionally, the method may include determining that a collision probability exceeds a threshold, the collision probability representing a probability of a collision between a subject Earth-orbiting satellite and a satellite causing one of the satellite trails for which orbital parameters are calculated, based on aforward propagation in time of the orbit represented by the orbital parameters, and in response to determining that the collision probability exceeds the threshold, outputting an alert to a control system of the subject Earth-orbiting satellite.

[0049] Optionally, the method may include, at a control system of a controlled Earth orbiting satellite, determining that a collision probability exceeds a threshold, the collision probability representing a probability of a collision between a subject Earth-orbiting satellite and a satellite causing one of the satellite trails for which orbital parameters are calculated, based on a forward propagation in time of the orbit represented by the orbital parameters, and in response to determining that the collision probability exceeds the threshold, outputting a control signal to the subject Earth orbiting satellite to modify the orbit thereof.

[0050] Embodiments may include apparatus comprising processor hardware and memory hardware, the memory hardware storing processing instructions which, when executed by the processor hardware, cause the processor hardware to perform a method including: obtaining image data representing a field of view of an optical system, the image data representing one or more satellite trails and a celestial background, and being mapped to a coordinate system in fixed spatial relation to the celestial background; executing differential image processing on the obtained image data to extract one or more satellite trails from the celestial background; executing matched filtering to calculate orbital parameters representing each of the one or more extracted satellite trails. The apparatus may further comprise the optical system.

[0051] Embodiments may include a computer program comprising processing instructions which, when executed by a computing apparatus comprising memory hardware and processor hardware, causes the computing apparatus to perform a method as set out above or elsewhere in the present disclosure. As a summary, embodiments provide one or more from among:

[0052] -a geometric framework for extracting satellite equatorial coordinates;

[0053] -a technique for mapping sensor pixels into the equatorial domain, including permitting the use of multiple sensors to form a synthesized image;

[0054] -an effective and computationally efficient technique for extracting transient events from the celestial background;

[0055] -an effective and computationally efficient matched-filter approach to calculate orbital parameters, demonstrated to be optimal under Gaussian noise;

[0056] -accurate orbital parameter calculations by modeling the optical system imaging the satellite trails and passing emulated trails through the model.Embodiments provide a framework for orbital estimation using optical sensors including wide field-of-view (Fo V) optical sensors, demonstrating their potential for space situational awareness. Embodiments decompose the orbital parameter estimation problem into three steps, significantly reducing the complexity compared with conventional techniques. By leveraging matched filtering techniques and differential image processing, embodiments successfully extract satellite trails from the celestial background and estimate their orbits. The efficacy of embodiments was validated through extensive sky simulations (set out in more detail below in the worked examples), showcasing the statistical performance of embodiments. Additionally, evidence of the efficacy of embodiments is set out in more detail below in a qualitative worked example based on real-world observations using the MASCARA observatory.

[0057] Embodiments may incorporate wide FoV sensors, which, when paired with the image processing techniques of the processes disclosed herein, detect and track satellites without prior knowledge of their orbits, offering a tool for enhancing space situational awareness.

[0058] Embodiments may include obtaining image data representing a field of view of an optical system, the image data representing one or more satellite trails and a celestial background, and being mapped to a coordinate system in fixed spatial relation to the celestial background; executing differential image processing on the obtained image data to extract one or more satellite trails from the celestial background; and executing optimization processing to calculate orbital parameters representing each of the one or more extracted satellite trails. The optimization processing may be executed by an algorithm such as particle swarm optimization.

[0059] Differential processing (see S102 explained below) and matched filtering (see S103 explained below) processes of embodiments enable satellite detection and tracking without needing prior satellite trajectory information.

[0060] Embodiments may use ground-based wide FoV optical sensors which reduce costs associated with satellite deployment and radar infrastructure.

[0061] Embodiments may use airborne or spaceborne wide FoV optical sensors which eliminate dependence on weather conditions and light pollution.

[0062] Embodiments may provide real-time processing, that is, near-instantaneous orbital data, for use in collision avoidance and SSA.

[0063] By leveraging wide FoV optical sensors, embodiments provide extensive sky coverage for efficient satellite detection.By utilising a differential processing technique (at S 102 explained below), embodiments extract satellite trails from the background without requiring prior knowledge, crucial for identifying unknown satellites and space objects.

[0064] By utilising matched filtering with iterative optimization (at S103 explained below), embodiments provide accurate real-time orbital estimation by minimizing discrepancies between observed and predicted satellite paths.

[0065] Embodiments are able to process image data from flexible sensor configurations including across multiple sites.

[0066] Embodiments are adaptable to various altitudes and can be configured for tracking across LEO, MEO, and GEO.

[0067] Embodiments may be scaled or adapted with additional sensors or processing modules for enhanced tracking accuracy.

[0068] Embodiments provide real-time orbital data without requiring prior knowledge, enhancing SSA for collision avoidance.

[0069] Ground-based, passive optical sensors are more affordable than satellite or radar systems. The wide FoV and image data processing techniques enable broad sky coverage and precise tracking.

[0070] Embodiments avoid the large infrastructure requirements of traditional radar-based systems.

[0071] Embodiments provide satellite operators with real-time tracking and collision avoidance data (see, for example, explanation of S104).

[0072] Embodiments offer a strategic, cost-effective tool for monitoring space objects without reliance on expensive satellite-based infrastructure.

[0073] Embodiments may assist research institutions in tracking unknown or newly launched satellites, contributing to broader research and SSA capabilities.

[0074] Embodiments may provide information to satellite operators and manufacturers seeking to optimise satellite deployment and avoid collisions.

[0075] Embodiments may enable real-time tracking and collision avoidance, minimizing new debris and supporting sustainable orbital environments.

[0076] Embodiments may enhance SSA without additional satellite launches, reducing light pollution and protecting dark skies for observation.

[0077] Embodiments may offer an affordable, ground-based alternative to space-based tracking, benefiting smaller organizations and government agencies.

[0078] Embodiments may enhance satellite safety and longevity, which is essential for the growth of telecom and Earth observation sectors.

[0079] Embodiments may provide real-time orbital data, improving decision-making and operational efficiency for satellite operators.Embodiments may provide real-time space monitoring, which may support national and international security.

[0080] Embodiments may reduce risks to essential satellite services like communications and navigation, benefiting society’s connectivity and emergency response.

[0081] LIST OF FIGURES

[0082] Embodiments will now be described, by way of example, with reference to the accompanying drawings, in which:

[0083] Figure 1 illustrates a process;

[0084] Figure 2 illustrates a process;

[0085] Figure 3 illustrates an apparatus;

[0086] Figure 4 illustrates a process;

[0087] Figure 5 illustrates Keplerian orbital parameters;

[0088] Figure 6 illustrates satellite orbit projections;

[0089] Figure 7 illustrates a simplified geometric layout of the Earth’s shadow;

[0090] Figure 8 illustrates an example of the ideal trails of four satellites with different orbital parameters; Figure 9 illustrates an example of experimental validation of the Poisson distribution;

[0091] Figure 10 illustrates an example of the noise distribution of a commercial DSLR camera sensor at room temperature compared with an ideal Gaussian function;

[0092] Figure 11 illustrates simulated night sky;

[0093] Figure 12 illustrates an RA / DEC projection of the celestial background;

[0094] Figure 13 illustrates an exemplary outcome of differential image processing;

[0095] Figure 14 illustrates example results of differential processing from the simulated sky scene in Figure 11;

[0096] Figure 15 (a) & (b) shows an example outcome of matched filtering processing;

[0097] Figure 16 (a) - (c) shows emulated sky orbital estimation;

[0098] Figure 17 illustrates a polar projection of five sensors showing a single exposure;

[0099] Figure 18 (a) & (b) illustrates a subset of matching results; and

[0100] Figure 19 illustrates apparatus.

[0101] Figures 1 to 4: Introduction

[0102] Figures 1 and 2 illustrate processes for transforming image data in which satellites are visible by their trails against a celestial background into determined values of orbital parameters of the satellites. Steps S 101 to S103 are common to both processes. The process of Figure 2 differs from the process of Figure 1 insofar as it includes an additional processing step S104 of controlling a movement of an Earthorbiting satellite in dependence upon the determined values of the orbital parameters. Furthermore,Figure 2 includes an optional mapping step SlOla, which may optionally be included in the process of Figure 1.

[0103] Figure 4 illustrates a process in schematic form and also including some of the data transferred from, to, and between, various of the processing steps. The process of Figure 4 is an example of the process of Figure 1. In particular, S401a is an example of SlOla, S402 is an example of S102, and S103 is exemplified by steps S4031 to S4034.

[0104] Figure 3 is a schematic illustration of an apparatus for performing processes such as those illustrated in Figures 1, 2, and 4. The apparatus comprises a processing system 300 and optionally also an optical system 310. The processing system comprises a processor 302 and a memory 304, along a data input / output interface 306 for receiving image data from the optical system 310. The optical system 310 may be part of the apparatus for performing the process, or may be excluded from the said apparatus, with the apparatus being configured to receive image data from the optical system 310 via an optical system interface 306.

[0105] Apparatus for performing the process may include an optical system, in which case the obtaining the image data may comprise an imaging process. Apparatus for performing the process may comprise computing hardware such as processor hardware and memory hardware in the absence of an optical system, such that the obtaining at S101 may comprise receiving the image data from an optical system.

[0106] The processor hardware 302 may comprise a CPU or a network of interconnected CPUs. Alternatively or additionally, the processor hardware may comprise a GPU or a network of interconnected GPUs. Alternatively or additionally, the processor hardware may comprise a NPU (neural processing unit) or a network of interconnected NPUs. Alternatively or additionally, the processing hardware may include an FPGA (Field-Programmable Gate Array) or a network of interconnected FPGAs. The or each processing unit may include a cache for temporary storage of data during processing tasks, with a control unit exchanging data and processing instructions with the memory hardware 304. The memory hardware 304 stores image data generated by the optical system 310 for processing by the processor hardware 302. The memory hardware 304 stores processing instructions which, when executed by the processor hardware 302, causes the processor hardware to execute a method such as illustrated in Figure 1. The processing system 300 of Figure 3 may be a computing device 900 such as is illustrated in Figure 19.

[0107] Step S101: ObtainingStep S101 comprises obtaining image data in which one or more satellites are visible as transient events or transient features. That is, the one or more satellites, by virtue of their motion, leave a transient feature in the image data with respect to an image background such as the celestial background. Step S101 may comprise obtaining image data representing a field of view of an optical system, the image data representing one or more satellite trails and a celestial background, and being mapped to a coordinate system in fixed spatial relation to the celestial background.

[0108] Obtaining the image data may comprise, at an optical system, taking one or a series of exposures of a field-of-view of the optical system by reading data from optical sensor hardware. The obtained image data may be images obtained by exposing the optical sensor hardware to the field of view, via optics, and reading from the optical sensor a signal representing an intensity of light received at each pixel. Alternatively, the obtained image data may comprise images captured by event-based camera hardware, the images being a series of difference images representing differences between a series of two or more exposures, and the differential image processing comprises computing a background image from the obtained series of difference images by calculating a mean image of the obtained difference images.

[0109] The terms transient feature and transient event are used interchangeably herein, and, unless specified to the contrary, refer to satellite trails, namely a visual artefact in an image or sequence of images representing a movement path with respect to a background of a satellite during the exposure of the image or sequence of images. The transient features may be continuous or may be composed of discrete portions (one per exposure) following the same orbital path.

[0110] The image data may be obtained by an optical system 310 comprising an optics 314 and an optical sensor 312. The optical system 310 may be airborne, spaceborne, or ground-based. The optical system 310 may comprise a ground-based optical sensor 312 on an equatorial mount, a fixed mount, or an AltAz mount. The optical system 310 may comprise an array of cameras, each camera comprising an optical lens 314 and an optical sensor 312 onto which incident light is focussed by the optical lens 314. An optical lens is an example of optics 314, which may also include, for example, a shutter. The optical system 310 also comprises a controller to read data out from the optical sensor to output to the processing system 300, and to control exposures of the optical sensor 312 by control of the shutter. The optical system 310 may be a camera or an array of cameras. The optics 314 may be a lens, a series of lenses, or a telescope.

[0111] The optical sensor hardware may be, for example, a wide FoV optical sensor such as a plurality of optical sensors arranged in an array. The optical system 310 may be satellite-mounted or otherwise spaceborne. The optical system 310 may be ground-based. The optical system 310 may be mounted onan equatorial mount. The optical system 310 may comprise an array of cameras, each camera comprising an optical lens and an optical sensor onto which incident light is focussed by the optical lens. The optical system 310 may comprise a plurality of event-based or neuromorphic cameras. The optical system 310 may comprise a single camera comprising an optical lens and an optical sensor onto which light is focussed by the optical lens. The single camera may comprise a fish-eye lens.

[0112] Orbital estimation for a Keplerian orbit involves obtaining the six parameters p = [a, e, i, fl, a>, v0], as explained in the below section relating to orbital projection. Step S 101 is a step of obtaining image data, from which satellite trails are to be extracted from the celestial background in the downstream processing at S 102. The image data may be obtained by an optical system 310. The celestial background is continuously moving relative to the sensor’s field of view (FoV). The optical system 310 may be configured to compensate for Earth’s rotation relative to the celestial background using equatorial mounts. An equatorial mount rotates the optical sensor around an axis aligned parallel to the Earth’s polar axis. The rotation rate is set to 360° per sidereal day to match the Earth’s average angular velocity. The equatorial mount may be equipped with a more sophisticated mechanism to implement frequent misalignment corrections. The optical system 310 may comprise a wide-FoV sensor. The optical system may comprise one or more sensors providing a combined FoV of greater than 30degrees, greater than 60degrees, greater than 90degrees, greater than 120degrees, or greater than 135degrees. Tracking via equatorial mounts imposes a practical limit on the FoV range, since the sensor would eventually intersect the horizon during rotation if not periodically readjusted.

[0113] As an alternative to mechanical tracking, a sequence of short exposures (for example, between 2 and 10 seconds per exposure, between 3 and 9 seconds per exposure, between 4 and 8 seconds per exposure) are taken and then digitally processed to compensate for the apparent motion of the stars in the mapping step SlOla discussed in more detail below.

[0114] The optical system 310 may comprise multiple cameras in an array formation, typically using silicon-based complementary metal-oxide-semiconductor (CMOS) or charge-coupled device (CCD) sensors 312, attached to optics 314 for focussing photons from a field of view onto the sensor 312. The optical sensor 312 may comprise a conventional RGB sensor in the visible spectrum, or infrared sensor (including near infrared, and thermal infrared), or neuromorphic sensor, or event-based sensor, or monochromic sensor.

[0115] Commercial off-the-shelf lens systems may be used which achieve FoVs in the range of 5°-190° with a single optical sensor 312. However, optical aberrations toward the edges of the image need to be corrected for. These aberrations include chromatic aberration due to increased optical dispersion forlarge-angle rays and coma aberration caused by the curved field projected onto typical flat sensors, which becomes more challenging to compensate for in large FoV systems.

[0116] As an alternative, the optical system 310 may comprise a camera array with multiple cameras and then apply digital signal processing to combine the images into a synthesized image. An example of such a system is MASCARA (Multi-site All-Sky CAmeRA), which comprises five cameras, providing an approximate total FoV of 135°. The optical system 310 may comprise a camera array with five cameras, which may be arranged to provide a combined FoV of greater than 90degrees, greater than lOOdegrees, greater than HOdegrees, greater than 120 degrees, greater than 130degrees, or greater than 140degrees. Given the known coordinates of the stars and with an accurate knowledge of time, in an initial calibration phase the pixels of the optical sensor hardware (whether a single sensor or sensor array) are mapped from the pixel domain naturally reported by the camera into the equatorial coordinate system, which may be referred to as the RA / DEC domain.

[0117] The optical sensor 312 may be sensitive to visible light (for example, light with wavelength in the range 380 to 700nm), infrared light (for example, light with wavelength in the range 700nm to 1mm), or nearinfrared light (for example, light with wavelength in the range 750nm to 1.4microns). Such image data may comprise readings only for pixels sensing a change in brightness or intensity of received light during an imaging period, with the readings quantifying the change.

[0118] S101a: Mapping

[0119] The extraction of transient features at S102 may be dependent upon the image data being mapped according to a particular coordinate system. That is, each pixel in the image data maps to a set of coordinates adhering to a particular coordinate system. In the processes of the present disclosure, including the processes of Figures 1 and 2, if the obtained image data at S 101 is expressed according to a coordinate system other than that required by the extract transient features step S102, then an intervening mapping step SlOla may be executed. That is, the mapping SlOla is a pre-processing step for the extracting transient features at S102. For the avoidance of doubt, the mapping step SlOla may be included in the process of Figure 1 if required.

[0120] The mapping step SlOla may be to map the image data, which may be referred to as optical observations, into the equatorial coordinate system, where the celestial background is considered static while satellites are highly dynamic. The obtained image data may be obtained by using astrometric tools and overlapping scenes collected from multiple sensors, for example.

[0121] In the specific example of S401a in Figure 4, the obtained image data is a sequence of raw image exposures in the pixel domain, which are processed by a coordinate transformation at S401a to obtaina sequence of projected image exposures. For example, the coordinate transformation is from the pixel domain into the equatorial coordinate system mapping each pixel to a 2-value coordinate comprising Right Ascension and Declination. S401a illustrates variable values and calibration parameters required for the coordinate transformation including UTC time, location, and camera calibration.

[0122] To achieve the mapping at S101a or S401a, the processing may comprise: (i) in an initial calibration, perform an astrometric calibration to find the Az / El coordinates for each pixel in every optical sensor 312 within the array (or single sensor), and (ii) for each time instance during the exposure or sequence of exposures, map the Az / El domain to the RA / DEC domain using coordinate transformation without needing to re-perform astrometric calibration for each subsequent image. An alternative approach is to apply astrometric calibration for each exposure and then map the pixels directly to RA / DEC coordinates; however, the first approach is more stable and requires less processing time than calibrating each image individually.

[0123] Since the images are now in the RA / DEC domain, the stars appear stationary, allowing for downstream processing at S102 onward.

[0124] An example of 100 images stacked in the RA / DEC domain is illustrated in Fig. 12. Fig. 12 illustrates RA / DEC projection of the celestial background obtained by averaging 100 exposures from the MASCARA sensor, each with an exposure time of 6.4 sidereal seconds (which is also an exemplary exposure time in embodiments). The stacking is achieved by taking the average across K = 100 images in the set λ, as follows:

[0125] sA

[0126] X P?. m) - ■■■■■;• V Afc (n,???). ( 17)

[0127]

[0128] Visibility of satellite trails is enhanced by specific processing techniques at SI 02, discussed in more detail below. In this example, the exposure time was set to 6.4 sidereal seconds, designed to ensure that the apparent motion of stars does not exceed a single pixel.

[0129] Orbital Projection

[0130] Embodiments may include mapping to the equatorial coordinate system from image data obtained by a ground-based optical system or by a satellite-based optical system. The present example relates to a ground-based optical system. Ideal satellite orbits are only affected by the homogeneous gravitational field and are thus independent of the Earth’s rotation. However, the apparent orbital observations made by a ground-based sensor are naturally influenced by the Earth’s rotational motion. In a spaceborne system, the projection follows the same steps except that the rate of background movement depends onthe attitude stability of the satellite rather than on the rotational speed of Earth. In typical applications onboard stable satellites, without an active attitude control system or with an attitude control system, the satellite exhibits very slow or controlled rotations relative to an inertial frame of reference. Accordingly, the exposure time per image can be extended beyond the values utilized on Earth's surface if required.

[0131] The coordinate system used to represent the satellite’s observed trail in the image data may be (i) the sensor’s local sky dome coordinate system, which uses elevation and azimuth angles referenced to the local horizon and the true north direction, or (ii) the equatorial coordinate system, typically employed for describing celestial objects independently of Earth’s rotation. Since the positions of background stars are well-documented and described in the equatorial system, astrometry methods can be used to precisely calibrate the sensor’s orientation and estimate any non-linear wide FoV lens effects. Additionally, as stars are relatively fixed compared to the rapid motion of satellites, subtracting the fixed background becomes more feasible and computationally efficient when using the equatorial coordinate system.

[0132] In the equatorial coordinate system, two angles are used to describe the positions of celestial objects: (i) the declination angle 8, which is the angle between the object and the equatorial plane (positive toward the north pole), and (ii) the right ascension a, which is the angle between the object’s spherical projection on the equatorial plane and the x-axis, pointing toward the vernal equinox (positive toward the east). Compared to satellites, such celestial objects are extremely distant and can thus be modeled as point sources located on a fixed sphere, known as the celestial sphere, which is an arbitrarily large sphere centred on Earth and aligned with its rotational axis, as shown in Fig. 6. The present example demonstrates a satellite coordinate propagation within the processes of Figures 1, 2, and 4, and shows how satellite trail observations are projected onto the celestial sphere by appropriate mapping, where required.

[0133] The following formulas provide a theoretical background for the necessary projections onto the celestial sphere and downstream processing. To construct an orbit, begin with an ideal equatorial ellipse, referred to as the base orbit, with the following Cartesian coordinates in the Earth-Centred Inertial frame.

[0134]

[0135] where v is the true anomaly angle and p is the instantaneous distance from the center, both of which are functions of time. The sine and cosine of the v were substituted in step (a) based on the known relation with the eccentric anomaly E. The distance p in ( 1) is calculated from the eccentric anomaly as follows: p = a(l — e cosE). The eccentric anomaly E is, in turn, calculated based on the mean anomaly M. The relation between the two is also well-known and provided as follows:

[0136] E

[0137]

[0138] = where M = f(E) = E — e sinE (2).

[0139] The mean anomaly M represents the orbital phase which linearly increases with time as follows:

[0140] M = (2π / T)t + M0(3)

[0141] where the time t is time variable, referenced to the satellite passing through the perigee, which thus coincides with the initial mean anomaly Mo, i.e. at t = 0. T is the orbital period, given by:

[0142] T = 2πa3 / 2 / √μ (4)

[0143] To calculate the initial mean anomaly Mofirst the initial eccentric anomaly Eois obtained by solving the following equation based on the provided initial true anomaly vo as follows:

[0144] ( tz... - R, 4RRR 1

[0145] R, ™ solution | an -- - - ™ - - -;— |, ( M

[0146]

[0147] [ 2 i - 4 R, |

[0148] where A = e / (l + Vl - e2) is a constant, accordingly Mo= Eo- e sin Eo. The based orbit is then subjected to three rotations respectively, (i) around the z-axis with an angle equal to the argument of periapsis o, (ii) around the x-axis with an angle i, and (iii) around the z-axis again with an angle Q, expressed as follows:

[0149]

[0150] where [x, y, z]Trepresents the Cartesian coordinates of the orbit as function of time t. For brevity, the propagation process can be abbreviated by a propagation function g(.), as follows:

[0151]

[0152] RR ™ O, (7) where s(t) = \x(t),y(t), z(t) is the state vector of the satellite, and p = [a, e, i, Ω, ω, v0] is the Keplerian parameter vector. The present example is predicated upon a simplified model of the Earth’s gravitational field. Processes may utilise a higher accuracy propagation function g(p, t) leveraging more complex methods such as the Simplified General Perturbation algorithm (SGP-4), which accounts for the Earth’s non-uniform gravitational field, atmospheric drag, solar radiation pressure, and the gravitational influence of the sun and the moon.

[0153] During the observation period, the ground sensor moves, in the ECI frame, due to Earth’s rotation, as indicated in Fig. 6, which shows the trail left by the sensor. Consequently, the apparent satellite trail, when projected on the celestial sphere, is no longer elliptical but is deformed according to the projectionseen from the sensor’s perspective. It is noted that a satellite-based optical sensor also moves and rotates during an exposure, and so by substituting the satellite’s movement’s for the Earth’s, the equivalent process for a satellite-based optical system is provided.Fig. 6 illustrates that satellite orbital projection on the celestial sphere does not depend on the eccentricity of the orbit.

[0154] In the ground-based example, let the sensor’s position vector be denoted as r(t). Then, the apparent position vector of the satellite is given by:

[0155]

[0156] Accordingly, the apparent equatorial coordinates of the satellite is found by converting the Cartesian vector into spherical angles as follows:

[0157]

[0158] atan2 is the two-argument arctangent function. Obtaining the ECI position of the observer, r(t) = [x0,y0,z0]T. is a process that first involves converting the observer’s geodetic coordinates, longitude 0. latitudeo, and height h0. into an Earth-Centered Earth-Fixed (ECEF) coordinate system. This is then converted into the ECI frame while knowing the Coordinated Universal Time (UTC). The conversion can be succinctly described as:

[0159] r

[0160]

[0161] (1) ™ HTC( 0 ) > (IQ) where UTC(t) represents the UTC date and time at a given perigee referenced time instance t.

[0162] Satellite Visibility

[0163] The processes of Figures 1, 2, and 4, rely upon satellite trails being present, that is, evident, in the obtained image data. For passive optical observation of Earth-orbiting satellites, there is an implicit reliance on the illumination provided by the Sun. Therefore, when a satellite passes through the Earth’s main shadow region, the umbra, solar light is completely eclipsed (see Fig. 7 for a simplified illustration showing the eclipsed segment of a satellite orbit). In low Earth orbits (LEO), due to the vast distance between the Earth and the sun, the shadow appears almost cylindrical, and the penumbra (the region where the Earth’s shadow partially eclipses the sun) is nearly negligible. The mathematical equivalent of determining the eclipsed segment of an orbit, i.e., the segment that lies within the umbra, can be expressed as follows:IEH! '):::<teci ($(£)» UTC(tl);sfi) € Umbra(UTC((E,

[0164] (11)

[0165] where Umbra(UTC(t)) is the cone representing the umbra at a given UTC date and time which can be numerically obtained based on the position of the sun in the ECI frame. Furthermore, in a typical scenario where the proposed trail sensor is pointed toward the zenith, the additional conditions required to detect a trail are: (i) the satellite must be above the horizon, and (ii) it must be within the optical sensor’s field of view (FoV). This requires the satellite to be above a minimum elevation angle threshold, θmin> 0. Accordingly, the visible section of the trail that falls within the FoV can be modeled as follows:

[0166]

[0167] where 0(t) is the satellite’s elevation angle. Note that if the sensor is oriented toward the zenith, the threshold is given by θmin= max(90° - FoV / 2.0). Consequently, the potentially visible section of a satellite trail is the logical intersection of all conditions:

[0168]

[0169] An example of a satellite simulation with respect to a ground-based sensor is shown in Fig. 8 for four different orbits. Note that one of the orbits is further eclipsed by the Earth’s shadow. Additionally, the visible RA / DEC region of the celestial sphere changes during the observation period due to the Earth’s rotation which is also indicated in Fig. 8. Fig. 8 illustrates an example of the ideal trails of four satellites with different orbital parameters, illustrating both the visible section of the sky and the eclipsed section of the orbit. The minimum elevation angle is set to θmin= 0.

[0170] S102: Extract Transient Features

[0171] Following obtaining the image data at S101 (and mapping into the equatorial coordinate system at SlOla if required), at S102 the obtained image data mapped in the equatorial coordinate system is processed to extract transient features or events. At S102 an image processing technique is executed to extract transient events from the static background, equivalent to extracting satellite trails. S102 may comprise executing image processing on the obtained image data to extract one or more satellite trails from the celestial background. The image processing may be differential image processing.

[0172] The image processing required at S102 is dependent upon the camera hardware and the obtained images. In a first example, the obtained image data comprises images captured by event-based camera hardware, the images being a series of difference images representing brightness or light intensity differences over a series of imaging periods, and the differential image processing comprises computing a background image from the obtained series of difference images by calculating a mean image of the obtained difference images. In a second example, the obtained image data comprises a sequence of exposures byoptical sensor hardware of the optical system executing the image processing comprises: (i) computing a series of difference images comprising a difference image between each pair of adjacent exposures in the sequence of exposures; (ii) computing a background image from the series of difference images by calculating a mean image of the series of difference images. In either example, the series of difference images are processed prior to the computing the background image by normalizing the series of difference images and / or bounding the pixel values to a predetermined saturation level.

[0173] Processes may utilise differential image processing to extract the transient features (i.e. the satellite trails) at S102. Since LEO satellite passes are transient, they are expected to traverse multiple pixels during the total exposure (image data acquisition) time, while the stars remain stationary. Accordingly, for a given subset I composed of K number of exposures, the first step is to compute a series of difference images between adjacent exposures in the sequence, this step is only required in non-event-based cameras while event-based cameras can provide the events directly from the sensor’s hardware. The difference between subsequent exposures may be computed as follows:

[0174]

[0175] where the max function across all exposures is performed in order to hold the trail across the K exposures, the subscript k indicates that the holding is performed across the exposures index k dimension. The processing at S102 which is applicable to both non-event-based cameras and eventbased cameras may further comprise normalizing the difference images. The processing at S102 may further comprise suppressing the image background. For example, the background may be suppressed as follows:

[0176] _zi _ ( | (19)

[0177]

[0178] max A \1 / 5

[0179] where the difference image is normalized with respect to the maximum value across the entire tensor ( RC'}

[0180] elements, and represents the normalized and bounded background of the subset 1 obtained as follows:

[0181]

[0182] is the normalized mean of the exposure subset I, and Pp(. ) represents the pthpercentile value. The purpose of bounding the background in (20) to a certain saturation level is to suppress the locations ofbright stars and mitigate their scintillation effects. The value of p is configurable to produce sharper trail images. is the mean image of the subset I, as indicated in equation (17).

[0183] Accordingly, a sequence of L processed subsets, referred to as frames, forms a new tensor as shown in Fig. 13, where each frame is the result of processing K exposures. Increasing the number of exposures forming a frame has two conflicting effects: (i) on one hand, background estimation in (20) improves due to the larger number of samples; (ii) on the other hand, satellite trails become longer, which reduces the temporal resolution of the resulting tensor. For low Earth orbits, an example number of frames could be in the range 4 to 16, while higher orbits such as medium Earth orbit could have more frames, e.g. range 16 to 50, due to the longer orbital period. Fig. 13 illustrates an exemplary outcome of the differential image processing for 5 frames of 25 simulated exposures each, along with the ground truth orbit.

[0184] SI 03: Determine Orbital Parameters

[0185] Following the extraction of the transient events from the background, the extracted transient events are the input for an orbital parameter determination step S103. S103 may comprise executing matched filtering to calculate orbital parameters representing each of the one or more extracted satellite trails.

[0186] The extracted transient events may be represented in image data or in some other data format in which the relevant characteristics are preserved (for example the pixel values of the pixels included in the trail and the equatorial coordinate system coordinates of the pixels, along with metadata such as UTC time).

[0187] A specific example is illustrated at S4031 to S4034. The orbital parameter determination step may comprise orbit estimation through a matched-filtering process S4031, which iteratively proposes S4032 orbital parameters to converge to optimal orbital parameters S4033 that minimize the discrepancy between observations and predicted orbits.

[0188] SI 03: Orbital Geometry

[0189] Orbital propagation is the process of estimating the future position and velocity of a satellite based on past measurements. The orbital parameters calculated at S103 may comprise values for each of six Keplerian elements, or for a subset thereof. The process is influenced by various factors, including the non-uniform gravitational field of Earth, the gravitational effects of the Moon and the Sun, solar radiation pressure, and atmospheric drag. Accurate orbital propagation typically requires the use of numerical methods to solve complex differential equations, as the full range of perturbing forces makes analytical solutions inaccurate over long periods. However, for the determination of orbital parameters in real-time with application for short-term use such as collision avoidance, an approximation can bemade by assuming a uniform Earth gravitational field and neglecting perturbations such as lunar and solar influences. Under these simplified conditions, the satellite’s trajectory follows an idealized elliptical orbit, which can be described using six parameters known as the Keplerian elements. These elements are examples of orbital parameters that may be determined at SI 03 (and combined with time and positional information, which may comprise initial true anomaly, to predict future position and velocity), and are briefly explained below in the context of the present disclosure:

[0190] - Semi-major axis: Denoted as a, represents the longest radius of the elliptical orbit. For a perfectly circular orbit, a would represent the radius of the orbit.

[0191] - Eccentricity: Denoted as e, describes how much the orbit deviates from being a circle. For typical satellite applications, the eccentricity is constrained by e e [0, 1) with e very close to zero represents a perfect circle, and e = 1 corresponds to a parabolic trajectory. Most Earthorbiting communication satellites aim for circular or near circular orbits to maintain relatively constant altitudes.

[0192] - Inclination: Denoted as i, is the angle between the normal of the orbital plane and Earth’s rotational axis (the z-axis).

[0193] - Right ascension of the ascending node: Denoted as fl, this angle represents the rotation of the orbital plane around the z-axis.

[0194] - Argument of periapsis: Denoted as m, describes the rotation of the ellipse within the orbital plane itself.

[0195] - Initial true anomaly: While the first five parameters are sufficient to describe the Keplerian orbit, the true anomaly specifies the satellite’s position along its orbit at a given time, called the epoch. This parameter is cmcial for orbital propagators to predict future positions based on the initial position and time elapsed since observation. It is denoted as v0.

[0196] The Keplerian elements are illustrated in Fig. 5, using the Earth-Centred Inertial (ECI) frame based on the International Celestial Reference Frame (ICRF). The z-axis of this frame is aligned with Earth’s rotational axis, and the x-axis points toward the vernal equinox. This reference frame is independent of Earth’s rotation, making the stars appear effectively fixed for the purpose of the extraction of satellite trails such as at SI 02 and S402.

[0197] SI 03: Optical sensor model

[0198] Determining the orbital parameters by matched filtering may comprise estimating a set of orbital parameter values for the one or more satellite trails in the image data and measuring a discrepancy between the obtained image data and a set of emulated satellite trails generated by the estimated set of orbital parameter values, and iteratively modifying the set of orbital parameter values to calculate an optimised set of orbital parameter values which satisfy a predefined optimisation condition. Since thetechnique outlined above is predicated upon a comparison between obtained image data and emulated satellite trails, processes at S103 may include enhancing the comparison by passing the emulated satellite trails through a model representing an optical system via which the image data is generated.

[0199] When performing astronomical imaging, photons are collected by individual pixels in the optical sensor hardware, which may be a single sensor or an array of sensors, over a given period, called the exposure time. The exposure time may be a single continuous exposure, or may be a series of discrete exposures. Due to the inherent quantum nature of photons, the total number collected during the exposure is randomised, typically following a Poisson distribution. Through the photoelectric effect, incident photons release electrons in the sensor, which accumulate in each pixel well within the exposure time. This build-up of electric charge is then read out using an amplification circuit and converted into a digital signal via analog-to-digital converters. A monochromatic digital image is subsequently stored as matrix!(n, m), where n and m are the sensor’s pixel coordinates, and I represents the intensity at each pixel along with other sources of noise. Accordingly, if the intensity generated by the light falling on the sensor is denoted as I]ight, it may be expressed as follows:

[0200]

[0201] ~? / («•, m) m): E ~ Pois, (14) with r) (n, m) representing the efficiency of the pixel (n, m) in converting light into photoelectrons, Pois (Aught) is the Poisson distribution with mean Alight, and E denotes the collected optical energy from the object, which is proportional to the light flux density at pixel

[0202]

[0203] A key feature of the Poisson distribution is that its variance is equal to its mean. An example of this relationship is shown in Fig. 9, which depicts an actual stack of astronomical images captured with a DSLR camera sensor. Fig. 9 illustrates an example of experimental validation of the Poisson distribution including in (a) a stack of 41 astronomical images and in (b) variance of each pixel plotted against its mean value in arbitrary intensity units. The images were captured using a 203 mm narrow-filed telescope with a commercial DSLR camera.

[0204] In addition to the electrons generated by the incident light, there are additional electrons that are thermally agitated and independent of the incident light. These electrons cause what is known as thermal noise. Furthermore, the readout circuitry introduces additional noise. Therefore, a simplified model to represent the total intensity of a pixel can be expressed as follows:

[0205] $

[0206]

[0207] . h * T’iread ^bias-here, Ithermal and Iread represent the thermal and readout noises, respectively. The last term, Ib,as, is added by the readout circuitry to ensure that the output remains positive. Note that the image of the stars is convolved with the point spread function of the camera, denoted as h. The combined effect of noise and bias can be modeled as a Gaussian distribution, as shown in Fig. 10, which is based on an actual Canon450D sensor utilised in testing the processes disclosed herein. Fig. 10 illustrates an example of the noise distribution of a commercial DSLR camera sensor at room temperature compared with an ideal Gaussian function. For short-term exposures, the impact of atmospheric scintillation may also be considered and corrected for, as rapid turbulence in the Earth's atmosphere causes fluctuations in the apparent intensity of stars. However, in the model presented herein, only the effect of shot noise using the Poisson distribution is included, for simplicity. Processes may be configured with a model compensating for rapid turbulence in the Earth’s atmosphere, as a further factor influencing the observed intensity of a star.

[0208] SI 03: Exemplary Orbital Parameter Determination Workflow

[0209] A specific example of orbital parameter estimation (which may be referred to as determination, calculation, computation) is illustrated in Figure 4 at S4031 to S4034. Matched filtering may be used as an approach to estimate the orbital parameters. Executing matched filtering may comprise estimating a set of orbital parameter values for the one or more satellite trails in the image data and measuring a discrepancy between the obtained image data and a set of emulated satellite trails generated by the estimated set of orbital parameter values, and optimising the estimated set of orbital parameter values through iteration to satisfy a predefined optimisation condition.

[0210] Orbital parameter calculation based on satellite trail observation is an inverse problem where a set of ground truth orbital parameters p = [a, e, i, fl, <n,v0] produces a noisy optical observation A = {A;, V I £ [0, L - 1]}. The processing at S103 estimates a vector p that minimizes the discrepancy between the actual optical observation and the emulated trails A (that is, the trails that would be left by a satellite or satellites orbiting according to the estimated orbital parameters at the exposure time(s)) based on the orbital and optical models. The estimate orbital parameters are iteratively adjusted to reduce a cost function representing discrepancy between the extracted trails from the obtained image data and a set of hypothetical or emulated trails that would be observed at the optical system by a satellite or satellites orbiting according to the estimated orbital parameters.

[0211] One approach to optimising the estimation of orbital parameters is to minimize a cost function gcostC, •) that represents the discrepancy, as follows:

[0212] p = min ( A. A H ' (2D

[0213]

[0214] p r / . I

[0215] Since A is a noisy measurement, the matched filter concept is implemented in S 102 since it is an optimal filter for maximizing the signal-to-noise ratio. A higher matched filter output corresponds to better matching, i.e. less discrepancy between the observation and emulation. The processor is configured to reduce the cost function. The cost function, gCOst can be described as follows:

[0216]

[0217] where ° denotes the element-wise multiplication of two matrices (Hadamard product). Executing the calculation of the cost function in (21) is computationally expensive, and alternatives to an exhaustive brute-force search may be more computationally efficient. For example, the iterative matched filtering S4031 and proposing orbital parameters S4032 may be processed in accordance with one or more optimization methods, such as the particle swarm optimizer.

[0218] Testing performance using the particle swarm optimiser has yielded favourable results, as set out below. An example of the estimation outcome is shown in Fig. 14, demonstrating the overlay of the same orbits presented in Fig. 11. Note that when multiple orbits exist in the same set of exposures an iterative detection process can be employed as follows: (i) apply the optimization S4031-S4032 in equation (21) to obtain p^. where q is the orbit index, being S4033 a single orbit estimation; (ii) generate a simulated trail of satellite q and at S4034 suppress it from the measurement ofA; (iii) repeat the iterative optimisation process S4031-S4032 for individual trails until a certain number of orbits Q is obtained, or until the cost function persists above a given threshold. Fig. 11 illustrates simulated night sky using a reduced Gaia star database, incorporating Poisson arrival effects, thermal noise, Gaussian point spread function, and projection non-linearity. The example also shows a 30-second trail of four satellites with an apparent magnitude of 6. Fig. 14 shows example results of differential processing from the simulated sky scene in Fig. 11, showing the detection of the four orbits.

[0219] Embodiments iterate through a process comprising (i) estimating parameters, (ii) measuring a cost function quantifying discrepancy between a projected appearance or shape of an emulated satellite trail that an orbit having the estimated orbital parameters would have in the image data and the observed satellite trail in the image data, (iii) determining whether the cost function measurement, or some other aspect of (i) and (ii), satisfies an optimisation condition (and if so, consider the parameters to be optimised). In cases in which there are multiple satellite trails in the image data, the process may also include eliminating the satellite trails from the image data for which the orbital parameters have been optimised, before progressing to optimising the orbital parameters for the next satellite trail in the image data.

[0220] An optimisation algorithm such as particle swarm optimisation may be used to estimate the parameters with the goal of minimising the cost function measurement. A search space for the set of orbital parameters is defined. The search space boundaries for each parameter may be based on theoretical extremes of the pertinent orbital parameter, or may be based on max / min or other bounds of the orbital parameter in a reference source such as a catalogue of Earth-orbiting satellites. Within that search space,particle positions may uniformly cover the search space. Smaller spacing between adjacent particles will increase the prospect of finding the optimum set of values (each particle corresponding to a set of orbital parameter values), but comes at a cost of computational processing overheads. The search space is n-dimensional, with n being the number of unknown orbital parameters. Optionally, n may be set to 6, for the six Keplerian elements. Optionally, one or more of the Keplerian elements (or other orbital parameters) may be set to a fixed value, which reduces the dimensionality of the search space accordingly.

[0221] SI 04: Collision Avoidance

[0222] The orbits represented by the orbital parameters calculated at S103 may be used as inputs to a collision avoidance system, that is, may be input to a controller of an Earth-orbiting satellite for use in collision avoidance procedures. S 104 may comprise determining that a collision probability exceeds a threshold, the collision probability representing a probability of a collision between a subject Earth-orbiting satellite and a satellite causing one of the satellite trails for which orbital parameters are calculated, based on a forward propagation in time of the orbit represented by the orbital parameters, and in response to determining that the collision probability exceeds the threshold, outputting an alert to a control system of the subject Earth-orbiting satellite. S104 may comprise determining that a collision probability exceeds a threshold, the collision probability representing a probability of a collision between a subject Earth-orbiting satellite and a satellite causing one of the satellite trails for which orbital parameters are calculated, based on a forward propagation in time of the orbit represented by the orbital parameters, and in response to determining that the collision probability exceeds the threshold, outputting a control signal to the subject Earth orbiting satellite to modify the orbit thereof.

[0223] An example of a collision avoidance method is detailed in Assman, Berger, & Grothkopp, The COLA Collision Avoidance Method, 2009, ESA, 5thEuropean Conference on Space Debris, vol 5, issue 1. An orbit of the controlled Earth-orbiting satellite is compared with the Q orbits represented by the orbital parameters by propagating the orbits forward in time for a time period of collision examination. A collision calculation is performed using the forward propagated orbits, including:

[0224] A close approach processing step to determine a point of time of minimum distance between the controlled Earth-orbiting satellite and each of the Q orbits (an example for how such an algorithm may function is presented in Salvatore Alfano (1994). Determining Satellite Close Approaches, Part II, Journal of the Astronautical Sciences, 42, 143;

[0225] A collision detection processing step determines, for each of the Q orbits and based on the determined minimum distances, whether the two objects (i.e. the controlled satellite and each of the Q satellites of the Q orbits) can collide at the time of closest approach by approximating a positional covariance ellipsoid representing position of the satellite at the pertinent time witha sphere (with a diameter equalling maximum ellipsoid expansion), and determining whether the spheres intersect.

[0226] If the spheres do intersect, a collision probability processing step calculates a collision probability is calculated. An exemplary technique is set out in Salvatore Alfano (2007). Review of Conjunction Probability Methods for Short-term Encounters, AAS-07-148, which is based upon evaluating the positional covariance ellipsoids.

[0227] If the collision probability exceeds a threshold, an alert is output. For example, the alert may be a message to a human operator of the controlled satellite. Alternatively or in addition to the alert, a modification to the orbital parameters of the controlled satellite may be determined that will reduce the collision probability to below the threshold. The procedure for determining the modification may be iterative to calculate the effect on the collision probability and optionally also to perform collision detection processing and collision probability processing where required for the new orbital parameters and the Q satellites, and to re-iterate where the modified parameters cause collision probability above the threshold for any of the Q satellites.

[0228] Evaluate Performance

[0229] To evaluate the performance of the proposed model, extensive simulations have been conducted using realistic star databases, followed by a statistical assessment of its accuracy. Additionally, qualitative validation has been performed using real-world data from the state-of-the- art wide-field observatory MASCARA (Multi-site All-Sky CAmeRA), originally designed for exoplanet detection.

[0230] In testing to validate the performance of the processes of Figures 1, 2, and 4, a high-fidelity model was used not only for the orbital dynamics and imaging but also for the celestial background. For creating the background star locations were sourced from the European Space Agency’s Gaia archive, based on the Gaia space observatory mission (2013-2025). Since the archive contains more than 1.7 billion stars, only the brightest 500,000 stars were extracted, and those with an apparent brightness magnitude higher than 9 (which corresponds to fainter stars) were further excluded. For computational efficiency, the simulated sky image was constructed directly in the RA / DEC domain by summing the intensity of all stars within an RA / DEC cell (n, m), as follows:

[0231] Aijghvbcrc) ~ y d) 4- Assu(t>., d)> (W

[0232]

[0233] Vi <s,: <>■ nud C m wherestar(a, 5) is the star intensity at given RA / DEC coordinates, calculated based on the apparent magnitude M( a, 5) as follows: Astar(<z, 8) = k x io~o,4M(“'4. where k is an arbitrary gain value. The intensity of the satellites, denoted asstar(a, 5) depends on multiple factors, including surface area, surface materials, illumination angles, and spectral band. This simulation allows the injection of satellite trails at any desired apparent magnitude where the apparent equatorial coordinates of the satellites areobtained as previously described in relation to orbital projection. An example of a 30-second exposure is shown in Fig. 11 after adding the background sky, the trails of four satellites, thermal noise, and bias, as described in equation (15). Note that the point spread function used in this example is a simple Gaussian-shaped filter to simulate the impact of the imaging system.

[0234] The processes of Figures 1, 2, and 4 have been evaluated in simulations. Evaluation results have been obtained from both (i) a quantitative emulation of the night sky with known ground truth of artificially injected trails and (ii) a qualitative observation based on real-life exposures.

[0235] Worked Example: Simulated Image Data

[0236] In a worked example of the processes, simulated image data with satellites having known orbital parameters was generated and used as an input to evaluate and illustrate the functionality of the downstream processing steps. In particular, a set of 16 known orbits was generated with a RAAN spacing of 10° and an inclination spacing of 5°. All other parameters were kept constant to focus on the performance of fl and i. The orbits are assumed circular, with e = 0, and an altitude of h = 1200 km. For repeatability, the initial simulation time vector 14 / 07 / 2024:00:00:00 was used, while the remaining parameters are listed in Table I. The exposure session consists of K = 120 exposures, equally spaced at 6-second intervals (shutter speed is the same as the exposure time as there are no gaps between exposures).

[0237] After generating the simulated image data comprising star background and satellite trails as explained above, the simulated exposures were grouped into frames, each composed of 12 exposures processed together using the processing method described above in relation to S 102. The matched filter technique of S103 was applied with iterative orbit extraction to estimate the 16 orbits. Given the circular orbit assumption, the orbital parameters eccentricity and argument of periapsis were known (e = 0, a> = 0) during the optimization process; all other Keplerian orbital parameters were assumed unknown and thus were calculated, i.e., [a, i, D, v0]. To obtain statistically viable outcomes, the estimation was repeated for 800 different orbits (in sets of 16). Each run iteratively generated new simulated image data comprising a set of 12 exposures and then the processing steps S102 and S103 (or equivalents from the Figure 4 process) executed to extract orbital parameters for the 16 orbits. An example outcome of one run is depicted in Fig. 15, showing the successful extraction of all the 16 orbits from one set. Fig. 15 shows an example outcome of a single run, including at (1) stacked sections of the emulator, and at (b) detected orbits sorted in order.

[0238] To determine the reliability of the orbital estimation, various methods can be used to detect outliers. One such method is the scaled Median Absolute Deviation (MAD), which is particularly robust againstoutliers. The cut-off deviation threshold was set to three times the scaled MAD, such that samples above this threshold were considered failures. Accordingly, the obtained retrieval success rate was 88.12%. For the successful retrievals, the resulting joint histogram of the RAAN and inclination errors, denoted as and 6j, respectively, is presented in Fig. 16(a). As is observable in the Figure, the covariance in the presented setup shows a strong correlation between

[0239]

[0240] and e;, with a calculated Pearson correlation coefficient of 0.87. Additional histograms of the errors in the semi-major axis a and the initial true anomaly v0are shown in Fig. 16-(b) and (c), respectively. Fig. 16(a) to (c). Figure 16 illustrates emulated sky orbital estimation including (a) two-dimensional histogram of the joint RAAN and inclination errors; (b) histogram of the semi-major axis error; (c) histogram of the initial true anomaly error.

[0241] In processing image data at S102 and S103, regardless of whether the image data is live image data from optical systems, or simulated image data in a testing and evaluation environment, the image data is subject to differential processing per S102 or S402 to extract the satellite trails from the celestial background. The extracted satellite trails may each be a portion of a composite image obtained by combining image data obtained over a series of exposures. The series of exposures may be a frame comprising n exposures, or may be a series of m frames each comprising n exposures (wherein n is, for example, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, or greater than 12; and m is, for example, 2, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 15, 20, 25, 30, or more than 30). The satellite trail is extracted by differential image processing based on computing the difference between adjacent exposures with a max function to hold the trail across K exposures in a subset and with a bounded background to obtain a normalised mean of the subset of exposures. Each subset is a frame, and the obtained image data comprises m frames.

[0242] Each satellite trail in the image sequence data in turn is subject to matched filtering and successive estimation in order to find a set of orbital parameters best describing the extracted satellite trail, including accommodating for noise and other distortions introduced by the optical system (and optionally also for the Earth atmosphere, in the ground-based example. A cost function is defined to measure the discrepancy between the expected satellite trail in the image sequence - based on estimated orbital parameters - and the actual extracted satellite trail. The estimated orbital parameters are iteratively refined to minimize this discrepancy using a stochastic optimization technique, such as particle swarm optimization, differential evolution, gradient descent, or an improved variant of particle swarm optimization. Certain orbital parameters may be given fixed values in order to reduce the variables and therefore the number of iterations required to reach a termination condition. For example, eccentricity and argument of periapsis may be set to 0. The termination condition may be a predefined number of iterations, may be a predefined confidence level, or may be a predefined stability conditionsuch as cease processing when the aggregate change in the orbital parameter estimations across y iterations is less than a threshold, where y may be, for example, 2, 3, 4, or 5.

[0243] Different factors that could impact the accuracy include the brightness of the satellite and the length of the observable trail, optical system aperture, and atmospheric conditions. However, a key feature of stochastic optimization algorithms, such as the particle swarm optimization used in the present worked example, is that repeated execution of the algorithm can result in improved estimation. This is evident from the overall median errors: RAAN |

[0244]

[0245] | = 0.0015°, inclination |£i | = 0.011°, semi-major axis |ea| = 4.81 km, and initial true anomaly | eVo| = 0.75°. However, such increased accuracy comes at the cost of repeated optimization runs. A termination condition may be set, which when satisfied, ends the processing of step S103. For example, the termination condition may be a minimum delta between orbital parameters generated by subsequent processing iterations or a set of y processing iterations.

[0246] Worked Example: Image Data from MASCARA observatory

[0247] The primary objective of the MASCARA observatory was to explore new exoplanets based on slight changes in star apparent magnitude during planet transits. It operates by monitoring a very large section of sky, approximately 135°, simultaneously. MASCARA consists of two stations, one located in each hemisphere, and is equipped with five cooled cameras (ML11002) fitted with Canon 24 mmf / 1.4 USM L II lenses.

[0248] The present worked example obtained, as image data, a dataset of images and selected 100 samples to process. The recorded exposures commenced at 00:23:19 UTC on 24 / 10 / 2022, with each exposure lasting 6.4 sidereal seconds. The first step in processing the raw images is to calibrate the orientation of the five cameras and obtain the distortion polynomial of the lenses. To perform this task, the worked example used the popular astrometry tool (as described in D. Lang, D. W. Hogg, K. Mierle, M. Blanton, and S. Roweis, “Astrometry.net: Blind astrometric calibration of arbitrary astronomical images, ’’The Astronomical Journal, vol. 139, no. 5, pp. 1782-1800, 2010.), which produces the World Coordinate System (WCS) information — a common framework in astronomy used to map the pixel coordinates of an image to real-world celestial coordinates (RA / DEC).

[0249] With the knowledge of the observatory location and capture time, the worked example further translated the RA / DEC coordinates to Azimuth and Elevation (Az / El) so that each pixel has a known transformation to the (Az / El) coordinates. An example of a single exposure is shown in Fig. 17 after combining the pixels from the individual cameras and projecting them into the Az / El coordinate system. Note the clearly visible satellite trails during this short exposure. The remaining 99 images were notpassed through the astrometry process; instead, the worked example used the obtained calibration mapping to translate each pixel to the corresponding RA / DEC, given the time of each exposure and the known Az / El of each pixel.

[0250] The five images are then combined in the RA / DEC domain to produce a continuous projection, with an example previously shown in Fig.12, which further stacks the 100 images to demonstrate the concept. After this, the process typically follows the same workflow as S102 to S103 discussed above in the context of the simulated image data worked example. Specifically, the calibration and mapping, and differential image processing S102 are identical to those used for the simulated sky worked example. For step S103 however, given that the experiment is conducted on real satellite orbits, the worked example extracted the historical orbital element database matching the experiment datel. Instead of blindly optimizing for the best fit, the worked example directly drew the proposed orbital elements from the historical database. For a qualitative comparison, the estimation results are presented in Fig. 18(a), showing the detected orbits overlaid on the overall differential image. In Fig. 18(b), the same orbits are depicted along with their satellite catalog numbers.

[0251] Figure 19 is a schematic illustration of a hardware arrangement of an apparatus for performing a satellite monitoring method. The apparatus is a computing device 900 but optionally includes an optical system 310. The methods, processes, protocols, and techniques for satellite monitoring and collision avoidance described herein may be performed by apparatus having an arrangement such as illustrated in Figure 19. Apparatus having processor hardware and memory hardware described in the present specification may include one or more devices 900 having an arrangement such as illustrated in Figure 19. A plurality of such devices may be interconnected over a network such as a Local Area Network or the internet. A cloud service including performing one or more of the methods, processes, protocols, and techniques described in the present specification may be performed by one or more devices 900 having an arrangement such as illustrated in Figure 19.

[0252] The computing apparatus of Figure 19 comprises a plurality of components interconnected by a bus connection. The bus connection is an exemplary form of data and / or power connection. Direct connections between components for transfer of power and / or data may be provided in addition or as alternative to the bus connection.

[0253] The computing apparatus 900 comprises memory hardware 991 and processing hardware 993. Further components are optional according to the implementation requirements, including a network interface 995, input devices 997, and a display unit 999. The display unit 999 and the processing hardware 993may cooperate to implement a graphical user interface. The display unit 999 may be a touchscreen display unit.

[0254] The computing apparatus 900 may be a server or may be a laptop or other form of personal computer having a touchscreen display unit.

[0255] The memory hardware 991 stores processing instructions for execution by the processing hardware 993. The memory hardware 991 may include volatile and / or non-volatile memory. The memory hardware 991 may store data pending processing by the processing hardware 993 and may store data resulting from processing by the processing hardware 993.

[0256] The processing hardware 993 comprises one or a plurality of interconnected and cooperative CPUs for processing data according to processing instructions stored by the memory hardware 991. The processing hardware 993 may include one or more FPGAs, CPUs, GPUs, or NPUs.

[0257] A computing apparatus may comprise one computing device 900 according to the hardware arrangement of Figure 19, or a plurality of such devices operating in cooperation with one another. For example, as interconnected servers in a client:server arrangement.

[0258] A network interface 995 provides an interface for transmitting and receiving data over a network. Connectivity to one or more networks is provided. For example, a local area network and / or the internet. Connectivity may be wired and / or wireless. The network interface 995 may function as the optical system interface 306 of Figure 3 described above.

[0259] Input devices 997 provide a mechanism to receive inputs from a user. For example, such devices may include one or more from among a mouse, a touchpad, a keyboard, an eye-gaze system, and a touch interface of a touchscreen. Inputs may be received over a network connection. For example, in the case of server computers, a user may connect to the server over a connection to another computing apparatus and provide inputs to the server using the input devices of the another computing apparatus.

[0260] A display unit 999 provides a mechanism to display data visually to a user. The display unit 999 may display user interfaces by which certain locations of the display unit become functional as buttons or other means allowing for interaction with data via an input mechanism such as a mouse. A server may connect to a display unit 999 over a network.TABLE 1

[0261] SYMBOLS AND PARAMETERS

[0262] 1ABLK I

[0263] SYMBOLS AKO PAR EiBRS

[0264]

[0265]

[0266]

[0267] Symb Parameter Value

[0268] 73 Earth average radius 0371 km / .t Earth standard gravitatbnal parameter 33Ox Wuu Satellite orbit semi ma|or axis 733 brs e OtO eceemrtcity 0 i Orbit melmation (45^,00®) Q Right aseansbn ef the ascending node (220350) at Argument of periapsis 0*

[0269]

[0270] o* Isstsal tme saosiOy -MO?< Sate ite altitude

[0271]

[0272] «• ™ A ■{■ 120 feus

[0273] Avg-fiiwe Earth tadi 03 1 B BfevaiOs aagfe

[0274] Ojisisj Minimum etevutma angte 10”

[0275] A(, ) Pfu eeted image n RA / DEC pixels

[0276] «x Right ssosnshn angle

[0277] e matmn angle

[0278] n-s Senxofs longitude 7(044’ Ity’W

[0279] ( t Sensor's h bsib 2501025’3 t<<;Sensor's height 240

[0280] R’ ToO number of ex. ates

Claims

CLAIMS1. A satellite monitoring method, comprising:obtaining image data representing a field of view of an optical system, the image data representing one or more satellite trails and a celestial background, and being mapped to a coordinate system in fixed spatial relation to the celestial background;executing differential image processing on the obtained image data to extract one or more satellite trails from the celestial background;executing matched filtering to calculate orbital parameters representing each of the one or more extracted satellite trails.

2. The method according to claim 1, whereinthe calculated orbital parameters comprise values for each of six Keplerian elements.

3. The method according to any of the preceding claims, whereinexecuting matched filtering comprises estimating a set of orbital parameter values for the one or more satellite trails in the image data and measuring a discrepancy between the obtained image data and an emulated satellite trail generated by the estimated set of orbital parameter values, and iteratively repeating the estimating and measuring to optimise the estimated set of orbital parameter values to satisfy a predefined optimisation condition.

4. The method according to claim 3, whereinthe iteratively estimating and measuring to optimise the estimated set of orbital parameter values is performed by an optimisation algorithm.

5. The method according to claim 4, whereinthe optimisation algorithm is a particle swarm optimisation algorithm.

6. The method according to any of claims 3 to 5, wherein the emulated satellite trail is passed through a model representing an optical system via which the image data is generated before the discrepancy is measured.

7. The method according to any of claims 3 to 6, wherein as each estimated set of orbital parameter values is optimised the corresponding orbit is eliminated from the obtained image data.

8. The method according to any of claims 3 to 7, wherein the proposed values for the Keplerian elements is a matrix of values comprising, for each of the one or more extracted satellite trails, a matrix element per orbital parameter.

9. The method according to any of the preceding claims, wherein there are multiple satellite trails in the image data, and the executing the matched filtering comprises:(i) estimating a set of orbital parameter values for the first satellite trail in the image data and measuring a discrepancy between the obtained image data and an emulated satellite trail generated by the estimated set of orbital parameter values, and iteratively repeating the estimating and measuring to optimise the estimated set of orbital parameter values to satisfy a predefined optimisation condition; (ii) simulate a satellite trail of the optimised set of orbital parameter values in the image data;(iii) supress the simulated trail from the image data;(iv) repeat steps (i), (ii), and (iii) for each remaining satellite trail.(v) output the optimised sets of orbital parameters and / or the simulated satellite trails.

10. The method according to any of the preceding claims, whereinthe obtained image data comprises images captured by event-based camera hardware, the images being a series of difference images representing brightness or light intensity differences over a series of imaging periods, and the differential image processing comprises computing a background image from the obtained series of difference images by calculating a mean image of the obtained difference images.

11. The method according to any of the preceding claims, whereinthe obtained image data comprises images of a sequence of exposures of optical sensor hardware by the optical system; andexecuting the differential image processing comprises:(i) computing a series of difference images comprising a difference image between each pair of adjacent exposures in the sequence of exposures;(ii) computing a background image from the series of difference images by calculating a mean image of the series of difference images.

12. The method according to 10 or 11, wherein the series of difference images are processed prior to the computing the background image by normalizing the series of difference images and / or bounding the pixel values to a predetermined saturation level.

13. The method according to any of the preceding claims, whereinobtaining the image data comprises, at the optical system, taking one or a series of exposures of a field-of-view of the optical system by reading data from optical sensor hardware.

14. The method according to any of the preceding claims, whereinthe optical system is satellite -mounted or otherwise spaceborne, or wherein the optical system is airborne.

15. The method according to any of the preceding claims, whereinthe optical system is ground-based.

16. The method according to claim 15, whereinthe optical system is mounted on an equatorial mount.

17. The method according to any of the preceding claims, whereinthe optical system comprises an array of cameras, each camera comprising an optical lens and an optical sensor onto which incident light is focussed by the optical lens.

18. The method according to any of the preceding claims, whereinthe optical system comprises a single camera comprising an optical lens and an optical sensor onto which light is focussed by the optical lens.

19. The method according to claim 18, whereinthe optical lens of the single camera is a fish-eye lens.

20. The method according to any of the preceding claims, wherein the optical system is configured to generate image data via one or more optical sensors sensitive to light in:the visible light range; orthe infrared light range; orthe near-infrared light range.

21. The method according to any of the preceding claims, further comprisingdetermining that a collision probability exceeds a threshold, the collision probability representing a probability of a collision between a subject Earth-orbiting satellite and a satellite causing one of the satellite trails for which orbital parameters are calculated, based on a forward propagation in time of the orbit represented by the orbital parameters, andin response to determining that the collision probability exceeds the threshold, outputting an alert to a control system of the subject Earth-orbiting satellite.

22. The method according to any of the preceding claims, further comprising,at a control system of a controlled Earth orbiting satellite, determining that a collision probability exceeds a threshold, the collision probability representing a probability of a collision between a subject Earth-orbiting satellite and a satellite causing one of the satellite trails for which orbital parameters are calculated, based on a forward propagation in time of the orbit represented by the orbital parameters, andin response to determining that the collision probability exceeds the threshold, outputting a control signal to the subject Earth orbiting satellite to modify the orbit thereof.

23. Apparatus comprising processor hardware and memory hardware, the memory hardware storing processing instructions which, when executed by the processor hardware, cause the processor hardware to perform a method including:obtaining image data representing a field of view of an optical system, the image data representing one or more satellite trails and a celestial background, and being mapped to a coordinate system in fixed spatial relation to the celestial background;executing differential image processing on the obtained image data to extract one or more satellite trails from the celestial background;executing matched filtering to calculate orbital parameters representing each of the one or more extracted satellite trails.

24. The apparatus according to claim 23, the memory hardware storing processing instructions which, when executed by the processor hardware, cause the processor hardware to perform a method according to any of claims 2 to 22.

25. The apparatus according to claim 23 or 24, further comprising the optical system.

26. A computer program comprising processing instructions which, when executed by a computing apparatus comprising memory hardware and processor hardware, causes the computing apparatus to perform a method according to any of claims 1 to 22.

27. A non-transitory computer-readable medium storing the computer program of claim 26.

Images