diff --git a/scripts/PrintSourceRaDec.py b/scripts/PrintSourceRaDec.py
index 96eec0f..14b7596 100755
--- a/scripts/PrintSourceRaDec.py
+++ b/scripts/PrintSourceRaDec.py
@@ -3,19 +3,18 @@
 Utility program to print source Coordinates in RA,DEC that astrometry.net found
 use this with .rdls (rdls.fits) FITS files after solving an image
 """
-from pathlib import Path
-from astropy.io import fits
+
 import argparse
 
+from astrometry_azel.io import get_sources
+
 
 p = argparse.ArgumentParser(description="show RA DEC in .rdls files after solving an image")
 p.add_argument("fn", help=".rdls file from astrometry.net")
 p = p.parse_args()
 
-fn = Path(p.fn).expanduser().resolve(strict=True)
-with fits.open(fn, "readonly") as f:
-    radec = f[1].data
+radec_sources = get_sources(p.fn)
 
-print(radec.shape[0], "sources found in", fn, "with ra,dec coordinates:")
+print(radec_sources.shape[0], "sources found in", p.fn, "with ra,dec coordinates:")
 
-print(radec)
+print(radec_sources)
diff --git a/src/astrometry_azel/io.py b/src/astrometry_azel/io.py
index d81f44d..226f480 100755
--- a/src/astrometry_azel/io.py
+++ b/src/astrometry_azel/io.py
@@ -29,6 +29,28 @@
     loadmat = None
 
 
+def get_sources(fn: Path) -> xarray.Dataset:
+    """
+    read source (star) coordinates from .rdls file
+
+    Parameters
+    ----------
+
+    fn: pathlib.Path
+        .rdls file computed by astrometry.net solve-field (astrometry_azel.doSolve())
+
+    Returns
+    -------
+
+    radec: xarray.Dataset
+        RA, Dec of sources
+    """
+
+    fn = Path(fn).expanduser().resolve(strict=True)
+    with fits.open(fn, "readonly") as f:
+        return f[1].data
+
+
 def rgb2grey(rgb_img):
     """
     rgb_img: ndarray