Create format for Shapelets

Tigger/Models/Formats/ASCII.py

@@ -139,6 +139,19 @@ def get_ang_field(name, units=ANGULAR_UNITS):
     def getval(num, scale=1):
         return None if (num is None or len(fields) <= num) else float(fields[num]) * scale
 
+    def getval_list(num, scale=1):
+        return None if (num is None or len(fields) <= num) else [float(x) * scale for x in fields[num].split(',')]
+
+    def triangular_index(val):
+        print("IN TRIANGULAR INDEX")
+        n = 1
+        x = n
+        print("n, x", n, x)
+        while x < val:
+            n += 1
+            x += n
+        return n
+
     # now process file line-by-line
     linenum = 0
     format_str = ''
@@ -202,9 +215,12 @@ def getval(num, scale=1):
                     polpa_field, polpa_scale = format.get('pol_pa_rad', None), 1
             # fields for extent parameters
             extent_fields = [get_ang_field(x, ANGULAR_UNITS) for x in ('emaj', 'emin', 'pa')]
+            shapelet_fields = [get_ang_field(x, ANGULAR_UNITS) for x in ('sbetal', 'sbetam')]
             # all three must be present, else ignore
             if any([x[0] is None for x in extent_fields]):
                 extent_fields = None
+            if any([x[0] is None for x in shapelet_fields]):
+                shapelet_fields = None
             # fields for reference freq and RM and SpI
             freq0_field = format.get('freq0', None)
             rm_field, rm_err_field = get_field('rm')
@@ -307,16 +323,36 @@ def getval(num, scale=1):
                     if any([x is not None for x in spi_err]):
                         spectrum.spi_err = spi_err
             # see if we have extent parameters
-            ex = ey = pa = 0
+            ex = ey = pa = beta_l = beta_m = shapelet_coeffs_l = shapelet_coeffs_m = 0
             if extent_fields:
                 ex, ey, pa = [(getval(x[0], x[1]) or 0) for x in extent_fields]
                 extent_errors = [getval(x[2], x[3]) for x in extent_fields]
+            elif shapelet_fields:
+                (scoeffs_field, _) = get_field("shapelet_coeffs")
+                beta_l, beta_m = [(getval(s[0]) or 0) for s in shapelet_fields]
+                shapelet_coeffs_list = [float(x) for x in fields[scoeffs_field].split(',')] #getval_list(scoeffs_field)
+
+                list_len = triangular_index(len(shapelet_coeffs_list))
+                shapelet_coeffs = [[0.0 for _ in range(list_len)] for _ in range(list_len)]
+                max_val = 1
+                i_ind = 0
+                j_ind = 0
+                for val in shapelet_coeffs_list:
+                    shapelet_coeffs[i_ind][j_ind] = val
+                    i_ind -= 1
+                    j_ind += 1
+                    if i_ind < 0:
+                        i_ind = max_val
+                        j_ind = 0
+                        max_val += 1
             # form up shape object
             if (ex or ey) and max(ex, ey) >= min_extent:
                 shape = ModelClasses.Gaussian(ex, ey, pa)
                 for ifield, field in enumerate(['ex', 'ey', 'pa']):
                     if extent_errors[ifield] is not None:
                         shape.setAttribute(field + "_err", extent_errors[ifield])
+            if (beta_l or beta_m or shapelet_coeffs_l or shapelet_coeffs_m):
+                shape = ModelClasses.Shapelet(beta_l, beta_m, shapelet_coeffs)
             else:
                 shape = None
             # get tags

Tigger/Models/ModelClasses.py

@@ -442,6 +442,21 @@ def strDescErr(self, delimiters=('"', "x", "@", "deg"), **kw):
         err[0] * DEG * 3600, delimiters[0], delimiters[1], err[1] * DEG * 3600, delimiters[0],
         delimiters[2], round(err[2] * DEG), delimiters[3])
 
+class Shapelet(Shape):
+    typecode = "Sha"
+
+    mandatory_attrs = ["sbetal", "sbetam", "shapelet_coeffs"]
+    optional_attrs = dict(ex_err=None, ey_err=None, pa_err=None)
+
+    def getShape(self):
+        return self.ex, self.ey, self.pa
+
+    def getShapeErr(self):
+        err = [getattr(self, a + '_err', None) for a in self.mandatory_attrs]
+        if all([a is None for a in err]):
+            return None
+        return tuple(err)
+
 
 class FITSImage(Shape):
     typecode = "FITS"