Multi-Centroid Extraction Method for High-Dynamic Star Sensors Based on Projection Distribution of Star Trail

Abstract

To improve the centroid extraction accuracy and efficiency of high-dynamic star sensors, this paper proposes a multi-centroid localization method based on the prior distribution of star trail projections. First, the mapping relationship between attitude information and star trails is constructed based on a geometric imaging model, and an endpoint centroid group extraction strategy is designed from the perspectives of time synchronization and computational complexity. Then, the endpoint position parameters are determined by fitting the star trail grayscale projection using a line spread function, and accurate centroid localization is achieved through principal axis analysis and inter-frame correlation. Finally, the effectiveness of the proposed method under different dynamic scenarios was tested using numerical simulations and semi-physical experiments. The experimental results show that when the three-axis angular velocity reaches 8°/s, the centroid extraction accuracy of the proposed method remains superior to 0.1 pixels, achieving an improvement of over 30% compared to existing methods and simultaneously doubling the attitude measurement frequency. This demonstrates the superiority of this method in high-dynamic attitude measurement tasks.

1. Introduction

To meet the increasingly complex demands for spatial information acquisition, dynamic imaging mode has become the mainstream trend in the development of remote sensing technology [1,2,3]. Unlike traditional static imaging modes, dynamic imaging mode can overcome the contradiction between the field of view and the focal length of optical payloads, achieving high-resolution wide-field optical imaging [4,5,6]. The core of dynamic imaging mode lies in high-precision attitude control under maneuvering conditions, a characteristic that places more stringent requirements on the satellite’s attitude measurement capabilities during high-dynamic operations [7,8,9].

Star sensors are currently the most precise devices for attitude measurement in the aerospace field [10,11,12]. By capturing star images, they accurately extract the centroid positions of stars, using star image recognition and attitude determination algorithms to provide high-precision attitude measurements for spacecraft [13,14]. Under static imaging conditions, the measurement accuracy of star sensors can reach the arcsecond level, perfectly meeting the attitude measurement needs in static imaging [15].

However, the attitude measurement frequency of star sensors is limited by exposure time, making it impossible to provide high-frequency attitude measurement results [16]. In addition, under dynamic conditions, star spots no longer exhibit Gaussian distributions but instead form complex trajectories, known as star trails, whose shapes vary with coordinates and velocity. This increases the uncertainty in centroid localization, thereby reducing the accuracy of attitude determination [17,18,19]. These challenges severely limit the application of star sensors in dynamic optical imaging.

To address these issues, many studies have focused on improving the centroid localization accuracy and measurement frequency of star sensors under dynamic conditions. The main methods for centroid localization under static conditions include the grayscale weight method [20,21] and the Gaussian fitting method [22,23]. To enhance the centroid localization accuracy of star sensors in dynamic conditions, researchers have analyzed and compensated for centroid localization errors under dynamic conditions [24,25].

As the dynamic imaging models of star sensors become increasingly sophisticated, some imaging models have been applied to the centroid extraction of high-dynamic star images [26]. Some researchers have developed motion models for specific scenarios to deconvolve and restore star point trajectories, aiming to improve centroid positioning results by increasing the signal-to-noise ratio [27,28]. Some scholars have attempted to directly fit star point trails using prior information from established models to obtain the centroid coordinates [29,30]. However, these methods rely on the accurate modeling of the mission scenario, which limits their applicability.

In addition to improving accuracy, researchers have attempted to break through the limitations of imaging frame rates on measurement frequency. By controlling exposure, star trails are imaged as several short segments, and centroids are processed independently for different segments to estimate the attitude at different times within a single image, thus increasing the update rate of star sensors [31,32,33]. However, the implementation of multi-exposure star sensors often relies on electronic shutters or image intensifiers, not only increasing hardware complexity but also depending on the accurate estimation of motion parameters.

These two types of methods cannot simultaneously improve both the measurement frequency and centroid localization accuracy, making it difficult to meet the requirements of high-dynamic attitude measurements. In essence, star trails are the cumulative energy tracks of star points moving on the image plane during the exposure period, recording the changes in centroid positions over time [17,34]. Therefore, if this prior information can be used as a constraint to recover multiple centroid coordinates from star trails, it may be possible to improve both the measurement accuracy and the frequency of dynamic star sensors.

This paper proposes a new method for extracting centroids from star trails, which aims to simultaneously improve the centroid localization accuracy and double the measurement frequency of high-dynamic star sensors by extracting two precise endpoint centroid coordinates from the star trails. To achieve this, this paper first models the moving star image based on geometric imaging models and line spread functions, analyzes the projection distribution characteristics of star trails, and develops a centroid extraction strategy from the perspectives of time synchronization and algorithmic complexity. Subsequently, this paper details the algorithmic process for extracting endpoint centroids from star trails. Finally, this paper demonstrates an improvement in centroid localization accuracy and the frequency of high-dynamic star sensors through experiments, proving the correctness and application prospects of the proposed method.

The structure of this paper is as follows: Section 1 introduces the research background and significance, Section 2 provides the theoretical basis for the method proposed in this paper, Section 3 presents the specific centroid extraction process, Section 4 conducts numerical simulations and semi-physical experiments, Section 5 discusses the results, and Section 6 draws conclusions.

2. Models

To achieve centroid extraction from the star trails, it is necessary to model the motion of star points on the image plane during attitude maneuvers. For this purpose, this section, based on the geometric imaging model of the star sensor, establishes the mapping relationship between attitude information and star trail, analyzes the time synchronization of trajectory end points, and derives the projection distribution characteristics of the trails grayscale values. This provides a theoretical foundation for the subsequent methods proposed.

2.1. The Star Trail Model

As shown in Figure 1, the coordinate relationships in the star sensor are illustrated. Here, $O_I{X}_I{Y}_I{Z}_I$ represents the inertial space coordinate system. The star sensor is simplified to a pinhole imaging model with a focal length of $f$. The pinhole is taken as the origin $O_C$, and the $Z_C$ is the optical axis of the star sensor, establishing the camera coordinate system $O_C{X}_C{Y}_C{Z}_C$. By moving the $Y_P{X}_P$ plane along the $Z_C$ by $-f$, the $O_P{X}_P{Y}_P$ plane is obtained, where the $Z_C$ axis intersects with $Y_P{X}_P$ at $O_P$, establishing the two-dimensional image plane coordinate system $O_P{X}_P{Y}_P$. To simplify the description of pixel coordinates in image processing, the image plane coordinate system is translated by $\left[-W/2,-H/2\right]$ to obtain the image coordinate system $OXY$, where $W$ and $H$ are the width and height of the sensor, respectively.

At time $t$, the unit observation vector of the i-th star within the field of view of the star sensor in the camera coordinate system is ${\mathit{w}}_i(t)$

$${\mathit{w}}_i(t)={\displaystyle \frac{1}{\sqrt{x_i{(t)}^2+y_i{(t)}^2+f^2}}}\left[\begin{array}{c}-x_i(t)\\ -y_i(t)\\ f\end{array}\right]$$

(1)

where $\left[x_i(t),y_i(t)\right]$ is the centroid coordinates on the image plane, which are the horizontal and vertical coordinates. The unit vector of the star in the inertial coordinate system is ${\mathit{r}}_i$, described as

$${\mathit{r}}_i=\left[\begin{array}{c}\cos {\delta }_i\cos {\alpha }_i\\ \cos {\delta }_i\sin {\alpha }_i\\ \sin {\delta }_i\end{array}\right]$$

(2)

where ${\alpha }_i$ and ${\delta }_i$ are the right ascension and declination of the navigation star in the inertial coordinate system. The observation vector ${\mathit{w}}_i(t)$ and the star vector ${\mathit{r}}_i$ satisfy the following relationship

$${\mathit{w}}_i(t)={\mathit{M}}_i(t){\mathit{r}}_i$$

(3)

where ${\mathit{M}}_i(t)$ is the attitude matrix at time $t$, representing the transformation from the inertial coordinate system to the camera coordinate system. The attitude matrix can be defined in terms of Euler angles

$$\mathit{M}(t)={\mathit{R}}_z\left(\alpha (t)-{\displaystyle \frac{\pi }{2}}\right){\mathit{R}}_x\left(\beta (t)+{\displaystyle \frac{\pi }{2}}\right){\mathit{R}}_z\left(\gamma (t)\right)$$

(4)

where $\alpha$, $\beta$, $\gamma$ are the attitude angles of the camera coordinate system in the inertial coordinate system, and ${\mathit{R}}_z$, ${\mathit{R}}_x$ are rotation matrices.

Let the total exposure time of the star sensor be $T$, and the observation vectors of star ${\mathit{r}}_i$ at times $t_0$ and $t_1=t_0+\Delta t$ have the following relationship:

$${\mathit{w}}_i(t_1)=\mathit{M}(t_1)\mathit{M}{(t_0)}^T{\mathit{w}}_i(t_0)={\mathit{T}}_{t_0}^{t_1}{\mathit{w}}_i(t_0)$$

(5)

where ${\mathit{T}}_{t_0}^{t_1}$ is the state transition matrix from time $t_0$ to $t_1$. When the time interval $\Delta t$ is sufficiently small, the transition matrix approximately satisfies the following expression:

$${\mathit{T}}_{t_0}^{t_1}=\left[\begin{array}{ccc}1& \Delta {\xi }_z& -\Delta {\xi }_y\\ -\Delta {\xi }_z& 1& \Delta {\xi }_x\\ \Delta {\xi }_y& -\Delta {\xi }_x& 1\end{array}\right]$$

(6)

where $\Delta {\xi }_x$, $\Delta {\xi }_y$, and $\Delta {\xi }_z$ are the changes in the attitude angles around the three axes during the time interval at time $\Delta t$

$$\left[\Delta {\xi }_x,\Delta {\xi }_y,\Delta {\xi }_z\right]=\left[{\omega }_x,{\omega }_y,{\omega }_z\right]\Delta t$$

(7)

where ${\omega }_x$, ${\omega }_y$ and ${\omega }_z$ are the angular rates around the X, Y, and Z axes of the camera coordinate system, respectively.

Let $\left[x_i\left(t_0\right),y_i\left(t_0\right)\right]$ be the position of the centroid in the image coordinate system at time $t_0$, and $\left[x_i\left(t_1\right),y_i\left(t_1\right)\right]$ be the position of the centroid in the image coordinate system at time $t_1$. Based on Equations (1), (5) and (6), the displacement of the centroid between these two time points is

$$\left[\begin{array}{c}{x}_i\left(t_1\right)\\ y_i\left(t_1\right)\end{array}\right]=\kappa \left[\begin{array}{c}{x}_i\left(t_0\right)+f\Delta {\xi }_y+y_i\left(t_0\right)\Delta {\xi }_z\\ y_i\left(t_0\right)-f\Delta {\xi }_x-x_i\left(t_0\right)\Delta {\xi }_z\end{array}\right]$$

(8)

where $\kappa$ is the coefficient

$$\kappa ={\displaystyle \frac{1}{\left(-x_i(t_0)\Delta {\xi }_y+y_i(t_0)\Delta {\xi }_x\right)/f+1}}$$

(9)

where $x_i$ and $y_i$ are on the order of millimeters, $\Delta {\xi }_x,\Delta {\xi }_y$ are on the order of milliradians, and the focal length $f$ is more than 10 mm, thus within time $\Delta t$ one can take $\kappa =1$.

An ideal star sensor optically images distant stars onto the image plane, which can be considered as a Gaussian-type spot with energy distributed on a two-dimensional plane centered at the centroid.

$$E_s\left(x,y\right)={\displaystyle \frac{\eta \Delta t}{2\pi {\rho }^2}}\exp \left[-{\displaystyle \frac{{\left(x-x_i\right)}^2+{\left(y-y_i\right)}^2}{2{\rho }^2}}\right]$$

(10)

where $\rho$ is the radius of the Gaussian spot, and $\eta$ is a proportionality coefficient related to the energy of stars, the optical system, and the detector performance. During attitude maneuvers, the trajectory of the star point is represented by the accumulation of the Gaussian spot over the exposure time.

$$E_d\left(x,y\right)={\displaystyle \frac{\eta }{2\pi {\rho }^2}}{\displaystyle {\int }_{t_0}^{t_0+T}\exp \left[-\frac{{\left(x-x_i\left(t\right)\right)}^2+{\left(y-y_i\left(t\right)\right)}^2}{2{\rho }^2}\right]}dt$$

(11)

2.2. Time Synchronization of Endpoint Centroid Groups

The star sensor begins exposure at time $t_{start}$ and ends exposure at time $t_{end}$. Define the set of all centroids on the image plane at any time $t$ during the exposure period as a two-dimensional point set $\mathbb{W}\left(t\right)$, referred to as the centroid group at timestamp $t$

$$\mathbb{W}\left(t\right)=\left\{\left[x_i(t),y_i(t)\right]|i=1,2,\dots ,N\right\}$$

(12)

where $N$ is the number of centroids within the field of view. This paper assumes that $N>4$ meets the minimum star count requirement for star map recognition, and no centroids enter or leave the field of view during the exposure period. Attitude determination is then based on solving the three-axis attitude angles $\left[\alpha (t),\beta (t),\gamma (t)\right]$ using the centroid group $\mathbb{W}\left(t\right)$.

When the attitude is stable, centroids do not move during the exposure, and it can be considered that there is only one centroid group, which corresponds to time $t_{mid}=\left(t_{end}-t_{start}\right)/2$. However, as shown in Figure 2, during attitude maneuvers, the centroid coordinates become a function of time, resulting in a different centroid group corresponding to each moment during the exposure. When extracting the centroid group, it is necessary to consider not only its spatial position but also to determine the specific moment at which the centroid group corresponds.

To address this, this paper considers two specific moments during the exposure process to eliminate the time uncertainty of the centroid groups. These are the starting centroid groups $\mathbb{W}\left(t_{start}\right)$ and ending centroid groups $\mathbb{W}\left(t_{end}\right)$ corresponding to the start and end of the exposure, which are located at the starting and ending points of each star point trail. Clearly, the start and end times of the exposure are controlled and well defined, while the coordinates of the centroids are retained at the two end points of the star trails.

2.3. Projection Distribution Characteristics of the Star Trail

Although Equation (11) provides a model for the star trail, the starting and ending centroid coordinates appear in the model as time-varying parameters, making it extremely difficult to directly solve. To simplify this solution process, this paper proposes using the projection distribution characteristics of the star trail as a prior, converting the two-dimensional time-dependent fitting into two one-dimensional fittings.

From Equation (8), it is known that the displacement of the centroid can be decomposed into components along the X-axis and Y-axis of the image coordinate system. Taking the X-axis direction as an example, the displacement consists of two parts: one is a constant displacement $f\Delta {\xi }_y$ caused by the angular velocity motion, which is proportional to the focal length of the star sensor; the other is a time-varying displacement $y_i\left(t_0\right)\Delta {\xi }_z$ caused by the motion around the optical axis, which is proportional to the Y-axis coordinate of the centroid.

Considering a practical imaging system, the focal length $f$ is typically 1~2 orders of magnitude larger than the image plane dimensions. Therefore, the constant displacement is the dominant part of the centroid displacement. This means that during the exposure period, the displacement speed of the star points along the X-axis and Y-axis can be considered constant.

To study the uniform motion process of star points on the image plane, we first start by discussing the motion of a one-dimensional point spread function. Let the energy distribution of a static star point be represented by the $f_{PSF}$.

$$f_{PSF}\left(x|x_1\right)={\displaystyle \frac{1}{\sqrt{2\pi }\rho }}\exp \left(-{\displaystyle \frac{{\left(x-x_1\right)}^2}{2{\rho }^2}}\right)$$

(13)

where $x_1$ represents the one-dimensional center coordinate of $f_{PSF}$, and $\rho$ represents the spot radius. The integral of the motion along the X-axis can be expressed as line spread function $f_{LSF}$

$$f_{LSF}\left(x|x_1,x_2\right)={\displaystyle {\int }_{x-x_1}^{x-x_2}{f}_{PSF}\left(u\right)\mathrm{d}u}={\displaystyle \frac{1}{2}}\left[\mathrm{erf}\left({\displaystyle \frac{x-x_1}{\sqrt{2}\rho }}\right)-\mathrm{erf}\left({\displaystyle \frac{x-x_2}{\sqrt{2}\rho }}\right)\right]$$

(14)

where $x_1$ and $x_2$ represent the central coordinates at the start and end of $f_{PSF}$ movement, respectively, known as the start point and end point of the $f_{LSF}$. The term $\mathrm{erf}\left(\cdot \right)$ refers to the error function, which is the integral of the Gaussian function

$$\mathrm{erf}\left(x\right)={\displaystyle \frac{2}{\sqrt{\pi }}}{\displaystyle {\int }_0^x\exp \left(-t^2\right)dt}$$

(15)

Figure 3 illustrates the process of forming LSFs through the integration of PSFs with different spot radii.

Furthermore, the one-dimensional uniform motion of a two-dimensional star point is discussed. Let the start point of the two-dimensional Gaussian star point be $\left(x_1,y_0\right)$ and the end point be $\left(x_2,y_0\right)$ within the exposure time $T$. The energy distribution based on $f_{LSF}$ is

$$\begin{array}{c}{E}_m(x,y)={\displaystyle {\int }_{x-x_1}^{x-x_2}{f}_{PSF}\left(x\right)}{f}_{PSF}\left(y\right)dx=f_{LSF}(x)f_{PSF}\left(y\right)\hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}={\displaystyle \frac{1}{2\sqrt{2\pi }\rho }}\left[\mathrm{erf}\left({\displaystyle \frac{x-x_1}{\sqrt{2}\rho }}\right)-\mathrm{erf}\left({\displaystyle \frac{x-x_2}{\sqrt{2}\rho }}\right)\right]\exp \left(-{\displaystyle \frac{{\left(y-y_0\right)}^2}{2{\rho }^2}}\right)\hfill \end{array}$$

(16)

The two-dimensional grayscale distribution of the star point is shown in Figure 4. From Equation (16), it can be seen that the grayscale distribution is independent in the coordinate axis directions. On the velocity direction cross-section (parallel to the XOZ plane, shown in orange), the grayscale distribution is represented by $f_{LSF}$, and on the cross-section perpendicular to the velocity direction (parallel to the ZOY plane, shown in blue), the energy grayscale becomes $f_{PSF}$, i.e.,

$$\begin{array}{l}{E}_{mx}\left(x,y_i\right)=f_{PSF}\left(y_i\right)f_{LSF}(x)\\ E_{my}\left(x_i,y\right)=f_{LSF}(x_i)f_{PSF}\left(y\right)\end{array}$$

(17)

For a specific row, $f_{PSF}\left(y_i\right)$ and $f_{LSF}(x_i)$ are coefficients. Considering the summation of the star trail grayscale along the coordinate direction means calculating the marginal projection distribution of the star trail along the coordinate direction.

$$\begin{array}{l}{E_{px}}^{\prime }(x)={\displaystyle \sum _i^h{E}_{mx}\left(x,y_i\right)}=f_{LSF}(x){\displaystyle \sum _i^h{f}_{PSF}\left(y_i\right)}\\ {E_{py}}^{\prime }\left(y\right)={\displaystyle \sum _i^w{E}_{my}\left(x_i,y\right)}=f_{PSF}\left(y\right){\displaystyle \sum _i^w{f}_{LSF}(x_i)}\end{array}$$

(18)

where $w$ and $h$ are the width and length of the effective pixels of the star trail, respectively. The accumulation process only involves the superposition of coefficients, preserving the structural properties of $f_{PSF}$ and $f_{LSF}$, meaning that the centroid information of start point and end point is retained in the projection distribution.

If the image point has a two-dimensional uniform motion relative to the image plane, let the start point of the image point be $\left(x_1,y_1\right)$ and the end point be $\left(x_2,y_2\right)$. As shown in Figure 5, the grayscale distribution of the star trail is no longer parallel to the coordinate axes, and the velocity of the star point on the image plane is decomposed into $v_x$ along the X-axis and $v_{\mathrm{y}}$ along the Y-axis.

The cumulative projections of the star trail along the two axis directions are the result of the two component motions, as analyzed in the Y direction in Equation (18)

$$\begin{array}{l}{E}_{px}(x)={\displaystyle \sum _i^h{E}_{mx}\left(x,y_i\right)}=f_{LSF}(x){\displaystyle \sum _i^h{f}_{PSF}\left(y_i\right)}\\ E_{py}\left(y\right)={\displaystyle \sum _i^w{E}_{my}\left(x_i,y\right)}=f_{LSF}(y){\displaystyle \sum _i^w{f}_{PSF}\left(x_i\right)}\end{array}$$

(19)

It can be seen that the projection distribution $E_{px}(x)$ of the star trail along the Y-axis and the projection distribution $E_{py}(y)$ along the X-axis still conform to the $f_{LSF}$.

In this paper, we only consider the situation where the projections of star trails on both the X-axis and Y-axis directions conform to a $f_{LSF}$. This implies that the star points must have a certain displacement in both directions. Therefore, it is necessary to impose constraints on the minimum displacement pixel $l_{\min }$ of star points in both directions [29], which is

$$\left\{\begin{array}{l}{\omega }_x{t}_{exposure}\ge 2\arcsin ({\displaystyle \frac{l_{\min }a{\cos }^2\varphi }{2f}})\\ {\omega }_y{t}_{exposure}\ge 2\arcsin ({\displaystyle \frac{l_{\min }a{\cos }^2\varphi }{2f}})\end{array}\right.$$

(20)

where $t_{exposure}$ represents the exposure time, $a$ is the pixel size, $\varphi$ is the half field-of-view angle. Assume $a$ is 5.5 μm, $\varphi$ is 8°, the focal length $f$ is 40 mm, and the exposure time $t_{exposure}$ is 100 ms. To achieve star trails greater than 20 pixels in both directions, the angular velocity on both axes needs to be greater than 1.6°/s.

Based on the above constraints, to solve the endpoint centroids of the star trail, one can fit its grayscale projection distribution using $f_{LSF}$ to obtain the four positional parameters that constitute the centroid coordinates. Although the final centroid coordinates cannot be obtained directly, the introduction of this prior significantly simplifies the calculation process compared to directly performing a time-dependent fitting on the two-dimensional distribution of the star trail. The determination of the endpoint centroid coordinates will be completed in conjunction with other features of the star trail in the following sections.

3. Methods

The modeling and analysis of star trails in the previous section indicated that the end points of star trails exhibit good time synchronization and grayscale distribution characteristics. This section will detail the multi-centroid extraction algorithm.

As shown in Figure 6, the high-dynamic star sensor begins exposure at time $t_{start}$ and ends exposure at time $t_{end}$, capturing a dynamic star image and establishing an image coordinate system $OXY$. The start coordinates and end coordinates of the grayscale distribution for each of the star trails is determined using the following three steps:

• Paramet er Fitting: The $f_{LSF}$ is used to fit the grayscale projections of the star trails along the X-axis, identifying 2 positional parameters ${\theta }_{x1}$ and ${\theta }_{x2}$ as the start or end positions of the star’s movement along the X-axis. The same process is repeated to the Y-axis to find ${\theta }_{y1}$ and ${\theta }_{y2}$.
• Coordinate Generation: Based on the angle $\phi$ between the principal axis of grayscale distribution and the X-axis, the 4 positional parameters from step 1 are combined into endpoint coordinates ${\left[x,y\right]}_{t1}$ and ${\left[x,y\right]}_{t2}$. Note that at this stage, it is still not determined which of $t_1$ and $t_2$ corresponds to $t_{start}$ or $t_{end}$.
• Timestamp Determination: By combining the results from the previous frame, the temporal order of the two endpoint coordinates ${\left[x,y\right]}_{t1}$ and ${\left[x,y\right]}_{t2}$ is determined, thus obtaining the ${\left[x,y\right]}_{start}$ and ${\left[x,y\right]}_{end}$ for the star trail.

Each ${\left[x,y\right]}_{start}$ and ${\left[x,y\right]}_{end}$ of all the star trails are then composed into centroid groups $\mathbb{W}\left(t_{start}\right)$ and $\mathbb{W}\left(t_{end}\right)$, respectively, completing the multi-centroid extraction for the dynamic star image.

3.1. Fitting Centroid Position Parameters of Star Trails

Let the dynamic star image $\mathcal{I}(u,v)$ have dimensions $U\times V$, and let $\mathcal{S}$ be the set of the star trail.

$$\mathcal{S}=\left\{{\mathcal{S}}_i|i=1,2,\dots ,N\right\}$$

(21)

where $i$ denotes the index of the star trail in the star image, and $N$ is the number of trails in the star image. The subset ${\mathcal{S}}_i$ represents the set of two-dimensional coordinates that belong to the connected region of a single star trail

$${\mathcal{S}}_i=\left\{\left(u,v\right)|\left(u,v\right) \mathrm{belongs} \mathrm{to} \mathrm{star} i\right\}$$

(22)

Calculate the projection distribution of a single star trail ${\mathcal{S}}_i$ in the row and column directions and normalize it.

$$\begin{array}{l}{P}_x(v)=Norm\left({\displaystyle \sum _{u=1}^U\mathcal{I}(u,v)}\right)\\ P_y(u)=Norm\left({\displaystyle \sum _{v=1}^V\mathcal{I}(u,v)}\right)\end{array}$$

(23)

where $P_x(v)$ and $P_y(u)$ are the projection distributions in the row and column directions, respectively.

To solve for the start and end points of the projection distributions, a fitting function $f\left(\mathcal{X};\theta \right)$ is designed based on one-dimensional $f_{LSF}$

$$f\left(\mathcal{X};\theta \right)={\theta }_A\left[\mathrm{erf}\left({\displaystyle \frac{\mathcal{X}-{\theta }_1}{\sqrt{2}{\theta }_{\rho }}}\right)-\mathrm{erf}\left({\displaystyle \frac{\mathcal{X}-{\theta }_2}{\sqrt{2}{\theta }_{\rho }}}\right)\right]$$

(24)

where $\theta =[{\theta }_A,{\theta }_1,{\theta }_2,{\theta }_{\rho }]$ are the parameters to be fitted, with ${\theta }_A$ being the amplitude of $f_{LSF}$, ${\theta }_1$ and ${\theta }_2$ being the start and end point coordinates of $f_{LSF}$, and ${\theta }_{\rho }$ being the spot radius. Parameters ${\theta }_1$ and ${\theta }_2$ will form the endpoint coordinates, which are of concern, while ${\theta }_A$ and ${\theta }_{\rho }$ are process variables intended to facilitate the convergence of the function and will not be of further interest to us subsequently.

Since it is necessary to perform two fits of the grayscale distribution projections along the X direction and the Y direction separately, ${\theta }_x=[{\theta }_{xA},{\theta }_{x1},{\theta }_{x2},{\theta }_{x\rho }]$ and ${\theta }_y=[{\theta }_{yA},{\theta }_{y1},{\theta }_{y2},{\theta }_{y\rho }]$ are used to represent the two sets of parameters to be fitted, respectively. Based on the mean squared error, the error functions in the row direction $erro{r}_x\left({\theta }_x\right)$ and the column direction $erro{r}_y\left({\theta }_y\right)$ are defined as follows:

$$\begin{array}{l}erro{r}_x\left({\theta }_x\right)={\displaystyle \frac{1}{V}}{\displaystyle \sum _{v=1}^V{\left(P_x(v)-f(v;{\theta }_x)\right)}^2}\\ erro{r}_y\left({\theta }_y\right)={\displaystyle \frac{1}{U}}{\displaystyle \sum _{u=1}^U{\left(P_y(u)-f(u;{\theta }_y)\right)}^2}\end{array}$$

(25)

Using nonlinear least squares fitting to find the parameter combination that minimizes the error functions, the fitting process can be expressed as

$$\begin{array}{l}{\theta }_x=\arg \underset{\theta }{\min } erro{r}_x({\theta }_x)\\ {\theta }_y=\arg \underset{\theta }{\min } erro{r}_x({\theta }_y)\\ s.t.\hspace{1em}\left\{\begin{array}{l}0.9\le {\theta }_1\le 1\hspace{1em},\hspace{1em}0.5\le {\theta }_4\le 5\hfill \\ \begin{array}{l}\min \left(v\right)<{\theta }_{x1}<{\theta }_{x2}<\max \left(v\right) \mathrm{in} {\mathcal{S}}_i\\ \min \left(u\right)<{\theta }_{y1}<{\theta }_{y2}<\max \left(u\right) \mathrm{in} {\mathcal{S}}_i\end{array}\hfill \end{array}\right.\end{array}$$

(26)

where the amplitude of $f_{LSF}$ is constrained to be within $\left[0.9,1\right]$, and the spot radius is assumed to vary between 0.5 and 5 pixels. The start and end points of the trail should not exceed the range of coordinates in the subset ${\mathcal{S}}_i$. To accelerate the convergence of the fitting process, the initial values of parameters ${\theta }_x$ and ${\theta }_y$, denoted as ${{\theta }_x}^{(0)}$ and ${{\theta }_y}^{(0)}$, respectively, are set as follows:

$$\begin{array}{l}{{\theta }_x}^{(0)}=[1,5,\min (v),\max (v)]\\ {{\theta }_y}^{(0)}=[1,5,\min (u),\max (u)]\end{array}$$

(27)

An example of implementation is shown in Figure 7. The fitting results, ${\theta }_{x1}$, ${\theta }_{x2}$, ${\theta }_{y1}$, and ${\theta }_{y2}$ are the 4 parameters that constitute the centroid coordinates.

3.2. Generating Centroid Coordinates

Due to the dimensionality reduction in the grayscale distribution description achieved through projection statistics, the fitting process can only ascertain the numerical values comprising the endpoint centroids, resulting in ambiguity in the coordinates of the end points. This is illustrated in Figure 8 right. It is necessary to further combine the major axis directions of the star trail grayscale for a more accurate determination.

In this paper, the principal axis direction of the star trail is analyzed based on eigenvectors, and the weighted centroid $\left(\overline{x},\overline{y}\right)$ of the pixels within the star trail ${\mathcal{S}}_i$ is calculated using the moment of inertia method.

$$\overline{x}={\displaystyle \frac{{\displaystyle \sum {\displaystyle \sum \mathcal{I}(u,v)}}\cdot v}{{\displaystyle \sum {\displaystyle \sum \mathcal{I}(u,v)}}}}\hspace{1em},\hspace{1em}\overline{y}={\displaystyle \frac{{\displaystyle \sum {\displaystyle \sum \mathcal{I}(u,v)}}\cdot u}{{\displaystyle \sum {\displaystyle \sum \mathcal{I}(u,v)}}}}$$

(28)

Calculate the covariance matrix $C$

$$\begin{array}{l}C=\left(\begin{array}{cc}{\sigma }_{xx}& {\sigma }_{xy}\\ {\sigma }_{xy}& {\sigma }_{yy}\end{array}\right)\\ {\sigma }_{xx}={\displaystyle \frac{1}{N}}{\displaystyle \sum {\left(v-\overline{x}\right)}^2}\\ {\sigma }_{yy}={\displaystyle \frac{1}{N}}{\displaystyle \sum {\left(u-\overline{y}\right)}^2}\\ {\sigma }_{xy}={\displaystyle \frac{1}{N}}{\displaystyle \sum \left(v-\overline{x}\right)\left(u-\overline{y}\right)}\end{array}$$

(29)

and solve for the eigenvalues and eigenvectors of the covariance matrix $C$. Select the eigenvector $a=\left(a_x,a_y\right)$ corresponding to the larger eigenvalue; then, the angle $\phi$ between the principal axis direction and the X-axis is given by

$$\phi =\arctan \left(a_y/a_x\right)$$

(30)

Based on the angle of the principal axis, the discussion is categorized as follows, as shown in Figure 8.

• If the principal axis direction is in the 1st and 3rd quadrants, meaning that the signs of $v_x$ and $v_y$ are the same, then $\left({\theta }_{x1},{\theta }_{y1}\right)$ and $\left({\theta }_{x2},{\theta }_{y2}\right)$ are taken as the end points of the star trail.
• If the principal axis direction is in the 2nd and 4th quadrants, meaning that the signs of $v_x$ and $v_y$ are opposite, then $\left({\theta }_{x2},{\theta }_{y1}\right)$ and $\left({\theta }_{x1},{\theta }_{y2}\right)$ are taken as the end points of the star trail, i.e.,

  $$\left\{\begin{array}{l}\left[{\theta }_{x1},{\theta }_{y1}\right]\left[{\theta }_{x2},{\theta }_{y2}\right], \phi \in \left.\left(-{\displaystyle \frac{\pi }{2}},0\right.\right]\cup \left.\left({\displaystyle \frac{\pi }{2}},\pi \right.\right]\\ \left[{\theta }_{x2},{\theta }_{y1}\right]\left[{\theta }_{x1},{\theta }_{y2}\right], \phi \in \left.\left(-{\displaystyle \frac{\pi }{2}},-\pi \right.\right]\cup \left.\left(0,{\displaystyle \frac{\pi }{2}}\right.\right]\end{array}\right.$$

  (31)

Apply the processing procedures from Section 3.1 and Section 3.2 to each star trail, and agree that the coordinates that correspond to ${\theta }_{y1}$ form the centroid group $\mathbb{W}\left(t_1\right)$, while the coordinates that correspond to ${\theta }_{y2}$ form the centroid group $\mathbb{W}\left(t_2\right)$. Up to this point, the construction of the two centroid groups is completed.

3.3. Determining Timestamps for Centroid Groups

The start and end exposure timestamps for a dynamic star image come from the onboard computer’s time system, typically with microsecond-level precision. These two timestamps correspond to two constructed centroid groups, but the temporal order of these centroid groups has not yet been determined. This paper determines the temporal order of the centroid groups by analyzing inter-frame associations, thereby achieving the matching of centroid groups with their corresponding timestamps.

In the continuous star images captured by a star sensor under dynamic conditions, the centroids move along approximately linear paths on the image plane. Therefore, for centroid i, if $t_0$ is taken as the reference time and $\left[x{(t_0)}_i,y{(t_0)}_i\right]$ as the reference position, the distance between the star point and the reference position at time $t$ is positively correlated with the time difference $t-t_0$. Based on this principle, the temporal proximity can be inferred by analyzing the distances on the image plane.

Let the previous centroid group with a known timestamp, as determined by the star sensor, be $\mathbb{W}\left(t_0\right)$ with the timestamp $t_0$. For the i-th centroid in $\mathbb{W}\left(t_0\right)$, take its coordinates in $\mathbb{W}\left(t_1\right)$ and $\mathbb{W}\left(t_2\right)$, and calculate the Euclidean distances $d{(t_1)}_i$ and $d{(t_2)}_i$ from its position at $t_0$ to its positions at $t_1$ and $t_2$, respectively. These distances represent the displacement over the time intervals $t_1-t_0$ and $t_2-t_0$

$$\begin{array}{l}d{(t_1)}_i={‖\left[x{(t_1)}_i,y{(t_1)}_i\right]-\left[x{(t_0)}_i,y{(t_0)}_i\right]‖}_2\\ d{(t_2)}_i={‖\left[x{(t_2)}_i,y{(t_2)}_i\right]-\left[x{(t_0)}_i,y{(t_0)}_i\right]‖}_2\end{array}$$

(32)

For the two centroid groups, calculate the average displacement $D\left(t_1\right)$ and $D\left(t_2\right)$ of all centroids relative to its coordinates at time $t_0$

$$\begin{array}{l}D\left(t_1\right)={\displaystyle \frac{1}{N}}{\displaystyle \sum _{i=1}^{N}d{(t_1)}_i}\\ D\left(t_2\right)={\displaystyle \frac{1}{N}}{\displaystyle \sum _{i=1}^{N}d{(t_2)}_i}\end{array}$$

(33)

The centroid group with the smaller average displacement is also closer in time to the time $t_0$, as shown in Figure 9.

$$\left\{\begin{array}{c}\mathbb{W}\left(t_{start}\right)=\mathbb{W}\left(t_2\right),\mathbb{W}\left(t_{end}\right)=\mathbb{W}\left(t_1\right),\hspace{1em}D\left(t_1\right)>D\left(t_2\right)\\ \mathbb{W}\left(t_{start}\right)=\mathbb{W}\left(t_1\right),\mathbb{W}\left(t_{start}\right)=\mathbb{W}\left(t_2\right),\hspace{1em}D\left(t_1\right)<D\left(t_2\right)\end{array}\right.$$

(34)

To determine the timestamp of the current dynamic star images centroid group, the $\mathbb{W}\left(t_{end}\right)$ from the previous frame can be used as $\mathbb{W}\left(t_0\right)$. For the first dynamic star image where multiple centroids are extracted, there is no preceding $\mathbb{W}\left(t_{end}\right)$ to serve as a reference. Therefore, the centroid extraction can be performed on the previous star images using a classical centroid extraction method, such as the weighted centroid method (Equation (28)), and the resulting centroid group $\mathbb{W}\left(t_{mid}\right)$ can be used as $\mathbb{W}\left(t_0\right)$. Since $\mathbb{W}\left(t_{mid}\right)$ does not participate in the attitude information calculation, it does not require particularly high centroid extraction accuracy.

4. Experiment and Results

This paper employs numerical simulations and semi-physical experiments to validate the feasibility of the proposed method. In the numerical simulation experiments, the centroid localization accuracy of the proposed method is compared with that of the classical centroid localization method across various high-dynamic scenarios. The performance of centroid localization under noisy conditions is also tested, and the effect of multiple centroid extraction on doubling the attitude frequency is demonstrated. For the semi-physical experiments, a setup for capturing dynamic star images was constructed using a turntable and a star simulator, through which the effectiveness of the proposed method in real images was verified.

4.1. Numerical Simulation Experiment

To verify the impact of the three-axis angular velocity on the centroid extraction accuracy, simulations of star trail are conducted. The simulation parameters are as follows: focal length $f=42 \mathrm{mm}$, pixel size $a=5.5 \mathsf{\mu }\mathrm{m}$, number of pixels 2048 × 2048, exposure time 100 ms, simulation step $\Delta t=0.1 \mathrm{ms}$, quantization bit depth 8 bit. The numerical simulation of the star trail is based on Equation (11), with Gaussian-distributed noise added

$$E\left(x,y\right)=E_d\left(x,y\right)+\mathcal{N}(\mu ,\sigma )$$

(35)

4.1.1. Impact of Three-Axis Angular Velocity on Accuracy

In the previous analysis of Section 2.3, the Z-axis angular velocity was neglected. To verify the impact of three-axis attitude angles on centroid extraction, a simulation experiment was conducted to test the extraction accuracy under angular velocities of ${\omega }_x={\omega }_y=3^{\circ }/\mathrm{s}$, $5^{\circ }/\mathrm{s}$, $8^{\circ }/\mathrm{s}$, and ${10}^{\circ }/\mathrm{s}$, with the Gaussian noise variance set to 0.3.

To validate the effectiveness of the algorithm across the entire image plane, experiments were repeated with a star trail at different positions on the image plane, as shown in Figure 10, which include the center of the image (a), the right edge of the image (b), and the lower left corner of the image (c).

The proposed method is compared with two other extraction approaches: the weighted centroid method [26] and a fitting method [30] based on a dynamic image model. The former represents a classical method, while the latter is an improved method specifically strengthened for star trails.

To evaluate the precision of the centroid extraction process, the Euclidean distance $\epsilon$ between the true centroid position and the computed centroid position within the image coordinate system is utilized as a metric

$$\epsilon =\sqrt{{\left[x_s(t_{ref})-x\left(t\right)\right]}^2+{\left[y_s(t_{ref})-y\left(t\right)\right]}^2}$$

(36)

where $x\left(t\right)$ and $y\left(t\right)$ denote the centroid coordinates obtained at time t using the centroid extraction method, while $x_s\left(t\right)$ and $y_s\left(t\right)$ represent the true values from the simulation.

However, since different centroid extraction methods aim to determine the centroid coordinates at different target moments, the true values corresponding to these specific moments $t_{ref}$ should be selected for comparison. For the weighted centroid method and the fitting method, which aim to determine the centroid coordinates at the midpoint of the exposure time [26,30], the reference time $t_{ref}$ should be set to the $t_{mid}=\left(t_{start}+t_{end}\right)/2$. In contrast, the proposed method targets the determination of the centroid positions at both the start and end moments of exposure. Therefore, the reference time $t_{ref}$ should be set to either $t_{start}$ or $t_{end}$, accordingly.

The subplot in Figure 10 illustrates the process of this experiment. According to Equation (4), the centroid coordinates are calculated with a step size of 0.1 ms. For an exposure duration of 100 ms, this processes 1000 independent centroid coordinates, as shown by the red dots in the subplot. The accumulation of point spread functions centered on these centroid coordinates forms a star trail. Among these centroid coordinates, the coordinate at $t_{mid}$ is selected as the true value for comparison with the weighted centroid method and fitting method, while the centroid coordinates at times $t_{start}$ and $t_{end}$ are chosen as the true values for the proposed method, as shown by the green squares.

The experimental results are illustrated in Figure 11, where the horizontal axis represents the angular velocity in degrees per second, and the vertical axis shows the centroid extraction error in pixels. Specifically, particular attention is given to marking the angular velocity ±8°/s and the centroid extraction error of 0.1 pixels by gray dotted line. Different methods are distinguished by distinct colors. The proposed method is represented by orange and green curves, the weighted centroid method is indicated by black, and the fitting method is shown in brown. Note that, since the proposed method can extract two centroids from each star trail, each plot contains two results for it.

From the experimental results, it can be observed that as ${\omega }_z$ increases, the errors of all three centroid localization methods increase. However, the error of the proposed method is less sensitive to changes in ${\omega }_x$ and ${\omega }_y$. In conditions with ${\omega }_x$ and ${\omega }_y$ both reaching up to 10°/s, the average centroid extraction accuracy of the proposed methods shows a significant improvement, with enhancements of 70% over the weighted centroid method and 31% over the fitting method. In dynamic scenarios where ${\omega }_x={\omega }_y\le {10}^{\circ }/\mathrm{s}$,$\left|{\omega }_z\right|\le 8^{\circ }/\mathrm{s}$, the centroid extraction accuracy of the proposed method is better than 0.1 pixels.

4.1.2. Impact of Noise on Accuracy

To evaluate the performance of the proposed method under noisy conditions, star trails were simulated with different levels of Gaussian noise. The image plane positions were randomized, and a total of 30 simulations were conducted. The three-axis angular velocities are set to specific values 8°/s. The noise variance is varied from 0.1 to 2.0.

Figure 12 presents a boxplot that displays the central tendency and distribution of the centroid extraction error for the proposed method under varying variances of Gaussian noise. The horizontal axis denotes the increasing variance of noise, while the vertical axis represents the centroid extraction error in pixels. On each box, the black marker indicates the mean error, and the lower and upper edges of the box correspond to the 25% and 75% percentiles of the data distribution, respectively. The whiskers extend to the furthest points that are not considered outliers, with outliers individually marked by hollow circles.

The results indicate that the centroid extraction error increases with the increase in Gaussian noise variance. However, under conditions where the three-axis angular velocities are as high as 8°/s and the Gaussian noise variance is 2.0, the centroid extraction error does not exceed 0.25 pixels.

4.1.3. High-Frame-Rate Attitude Determination Experiment

The proposed method can extract two centroid groups from a dynamic star image, and these two centroid groups can be independently used for attitude determination, thereby doubling the attitude measurement frequency of the star sensor. To verify the effectiveness of the algorithm in doubling the attitude measurement frequency in a sequence of star images, a numerical simulation experiment was conducted in this section.

The experiment simulates a star sensor that operates at an imaging frame rate of 8 Hz and has an exposure duration of 100 ms. Within the simulation period of 1000 ms, it acquires a total of eight images. Stars with magnitude Mv5.5 from the Tycho-2 catalog are used, with initial attitude angles [2, 5, 0] deg and three-axis angular velocities [−8, −8, 8] deg/s.

Figure 13 shows the 1st, 3rd, 6th, and 8th dynamic star images out of the 8 images, captured during the intervals of 0–100 ms, 250–350 ms, 625–725 ms, and 875–975 ms, respectively. The images have been inverted, where deeper colors indicate higher energy levels. It is evident that the high-dynamic attitude maneuvers cause the star points to appear as distinct trails.

The proposed method was employed to process these dynamic star images in order to obtain multiple centroid groups. In each star image of Figure 13, the green and orange points, respectively, represent the starting centroids and ending centroids of the star trails extracted by the proposed method. Additionally, these four subplots display the endpoint grayscale distribution details of the same star, with the true values of the centroid coordinates marked by “×”. All the green points and orange points in each image constitute the starting centroid group $\mathbb{W}\left(t_{start}\right)$ and the ending centroid group $\mathbb{W}\left(t_{end}\right)$. For example, in Figure 13c, the centroid group $\mathbb{W}\left(725 \mathrm{ms}\right)$, composed of all the orange points, corresponds to the centroid positions at the end of the exposure for that star image, which occurs at 725 ms.

The centroid groups are applied to star image recognition and attitude determination algorithms [13], and the results of attitude determination are compared with the attitude data at corresponding moments to validate the effectiveness of the centroid groups. Figure 14 illustrates the high-frame-rate attitude determination results. From top to bottom, it shows the attitude angle calculation results ${\omega }_x$, ${\omega }_y$, and ${\omega }_{\mathrm{z}}$, respectively. The horizontal axis in the figure represents the time of the simulation in seconds, while the left vertical axis indicates angular velocity in degrees, and the right vertical axis shows the angular error $\Delta \omega$ in arcseconds. The green straight lines represent the attitude angles with uniform change. Green squares denote the attitude angles calculated from the starting centroid group $\mathbb{W}\left(t_{start}\right)$, and orange squares denote those calculated from the ending centroid group $\mathbb{W}\left(t_{end}\right)$. Black asterisks “*” indicate the attitude errors at the current moment.

This experiment shows that the principal axis analysis in Section 3.2 and the timestamp analysis in Section 3.3 ensure the coordinate positions and temporal sequence of the multi-centroid groups. The generated centroid groups satisfy the algorithmic requirements for subsequent attitude measurement tasks, resulting in superior accuracy with attitude errors of less than 8 arcseconds in both the X and Y directions, and less than 20 arcseconds along the Z-axis. Thanks to the extraction of dual centroid groups from a single star image, the proposed method achieves a measurement rate of 16 Hz at an imaging frame rate of 8 Hz, thereby doubling the efficiency of attitude information measurement, which demonstrates its potential for application in dynamic conditions where attitude information changes rapidly.

4.2. Semi-Physical Experiment

To verify the effectiveness of the proposed algorithm on actual situation, a semi-physical experimental system was designed to acquire dynamic star images that contain real noise.

As shown in Figure 15, the semi-physical experimental system for dynamic star images consists of the system controller, the drivers, a three-axis turntable, a star sensor, and a star simulator. The star simulator is fixed on a multi-freedom adjustment platform, and the star sensor is mounted on the three-axis turntable. The coordinate system is defined as shown in the Figure 15, with the Z-axis coinciding with the turntable’s axis of rotation. Initially, the entrance pupil of the star sensor is aligned with the exit pupil of the star simulator.

At the start of the simulation, the system controller drives the star simulator to display a static star image in a specified direction. Simultaneously, the turntable is controlled to move at a set angular velocity through the driver, and the star sensor is controlled to capture images. After controller receives the images from the star sensor, it preprocesses them and uses the proposed method to calculate the centroids and determine the attitude. Based on the attitude determination results, the angular velocity is calculated, and the calculated results are compared with the control parameters of the turntable.

The specifications of the star sensor are as follows: focal length $f=41.8 \mathrm{mm}$, pixel size $a=5.5 \mathsf{\mu }\mathrm{m}$, number of pixels 2048 × 2048, exposure time 166 ms, and principal coordinates $[u_0,v_0]=\left[1035.23,1133.77\right]$.

The angular velocities for the X-axis and Y-axis of the turntable are set to 5°/s. According to Rodrigues’ rotation formula, the resultant angular velocity is 7.07°/s. Therefore, the rotation angle ${\psi }_0$ during the exposure time is

$${\psi }_0={7.07}^{\circ }/\mathrm{s}\times 0.166s={1.17}^{\circ }/\mathrm{s}$$

(37)

The experimental results are shown in Figure 16 where three consecutive dynamic star images were selected for the test. Based on the proposed method, two centroid groups were extracted from each dynamic star image. The green and orange points in each image represent the starting centroid group and the ending centroid group, respectively. The attitude determination algorithm was then used to calculate the corresponding three-axis attitude angles for each centroid group, and the rotation angles were solved using Rodrigues’ theorem.

Ultimately, the recovered rotation angles $\psi$ from the three dynamic star images were 1.16°/s, 1.28°/s, and 1.19°/s, differing from the actual control values ${\psi }_0$ by 0.01°/s, −0.11°/s, and −0.02°/s, respectively.

5. Discussion

This paper proposes a novel centroid extraction method that extracts two centroid groups from a single star image and validates the precision and noise resistance of the proposed algorithm through experiments.

In Section 4.1.1, we compare the performance of the proposed method with that of traditional centroid extraction methods under high-dynamic conditions where angular velocities are present along all three axes. The results show that the proposed method achieves higher centroid extraction accuracy in complex dynamic scenarios. This is because, under the influence of three-axis angular velocities, star trails form curves with small curvatures, leading to significant errors in classical centroid localization methods, which rely on the overall features of the star trails. In contrast, the proposed method focuses on the edge features of the star trail end points, making it less sensitive to overall deformations. However, since the star trail model in this paper assumes uniform linear motion across the image plane, the error increases with the increase in Z-axis angular velocity. Future work could reduce this error through more refined modeling.

In Section 4.1.2, we evaluate the performance of the algorithm under different Gaussian variances, and the results show that the proposed method still performs well. The proposed method is based on fitting the projection distribution, and the summation process of gray values along the X-axis (or Y-axis) during the projection calculation enhances the algorithm’s robustness to noise.

By extracting two independent centroid groups from dynamic star images, the star sensor can perform two attitude measurements in a single imaging process. With the imaging frame rate remaining unchanged, the proposed method achieves a doubling of the measurement frame rate. This has significant potential for applications in high-dynamic scenarios where attitude information changes rapidly, as further demonstrated by the experiments.

The experiments in Section 4.1.3, which involve accurate multi-centroid extraction from consecutive star images, achieve a doubling of the attitude measurement frequency, validating the effectiveness of the proposed algorithmic procedure. The semi-physical simulation experiments in Section 4.2, based on accurate centroid extraction from real star images, complete the angular velocity estimation task for a single star image, expanding the application scenarios of star sensors in high-dynamic tasks.

6. Conclusions

To improve the centroid extraction accuracy and efficiency of high-dynamic star sensors, this paper proposes a multi-centroid localization method based on the prior distribution of star trail projections. This method successfully achieves high-precision extraction of the endpoint centroid coordinates of dynamic star trails by establishing an accurate mapping relationship between attitude information and star trails, thereby enhancing the centroid extraction accuracy and frequency of high-dynamic star sensors. Experimental results show that when the three-axis angular velocity reaches 8°/s, the centroid extraction accuracy of this method remains better than 0.1 pixels, and the attitude measurement frequency is doubled compared to the frame rate, demonstrating the superiority of this method in high-dynamic attitude measurement tasks.