## 4 Discussion  
We have demonstrated that by utilizing advanced image deconvolution techniques, specifically the RL algorithm tailored to the _NuSTAR_ PSF, we can significantly enhance the spatial resolution of _NuSTAR_ X-ray images. By choosing the correct model of the telescope's PSF and incorporating it into the deconvolution procedure, it is possible to resolve closely spaced X-ray sources that are indistinguishable from standard imaging techniques. Our  
Figure 4: Pulse phase-resolved deconvolution of simulated images. The pulse phase is associated with the spin period of the pulsating source (right). (**Top**) Deconvolved images of the two sources at different phases of the pulsating (right) source. The left point source remains at a fixed intensity irrespective of the spin phases. (**Bottom**) Pulse profile showing the variation of intensity with spin phase for both the static (blue) and the pulsating (red) source (see Section 3.3 for details). application of the RL iterative method to the case of the HMXB pulsars SXP 15.3 and SXP 305 (right panel in Figure 3) illustrates this approach. Observations of these nearly overlapping sources with _NuSTAR_ have been challenged by the instrument's broad and off-axis-dependent PSF, leading to blended images in which pulsations from one or the other source are detected at different epochs without the two sources ever being spatially resolved. By employing a refined PSF model accounting for energy dependence and off-axis variations inherent in the telescope's optics, we increased the effectiveness of the RL deconvolution algorithm. The ability to resolve such X-ray sources has significant implications for astrophysical research. For HMXBs like SXP 15.3 and SXP 305, accurately determining their individual fluxes, spectral properties, and temporal behaviors is essential for understanding their accretion processes, magnetic fields, and evolutionary states [29]. In crowded pulsar fields (as in the SMC), where multiple X-ray sources may be located in close angular proximity, enhanced resolution facilitates a more precise characterization of the source population and their individual contributions to the X-ray luminosity of the entire region. We also demonstrated the technique of phase-resolved RL deconvolution in simulated images of a pulsating and a non-pulsating source in close proximity (10''). This technique failed to yield clean images of the SXP 15.3/305 pair when the data were divided and stacked in four phase bins. The division resulted in decreased photon counts in each bin, and the S/N ratio also decreased sufficiently for spurious artifacts to appear due to the enhanced noise. Nevertheless, the proof of concept offered by the simulated data (Figure 4) demonstrates that in cases where the pulse period is known, phase-resolved deconvolution can detect the variations in the brightness of nearly overlapping sources in the sky. It is also important to note that the same technique can be extended to energy-resolved deconvolution in spectroscopic observations of high S/N ratio, which will enable different emission components (and mechanisms) to be spatially and temporally disentangled. Furthermore, the improved detection capabilities enabled by refined PSF modeling and deconvolution can also support studies of transient phenomena. Detecting outbursts, pulsations, or other temporal variations from faint sources requires not only temporal resolution but also the ability to spatially isolate the sources in order to avoid contamination from nearby objects. By enhancing the spatial resolution, we improve the S/N ratio for individual sources, enabling the possibility of more sensitive timing analyses. Despite the advancements in this field of research, certain challenges and limitations remain: (a) The RL deconvolution procedure is sensitive to inaccuracies in the adopted PSF model; any mismatch between the modeled PSF and the actual instrumental response can introduce artifacts or spurious features in the cleaned image.