Python QTable.add_index方法代码示例

# 需要导入模块: from astropy.table import QTable [as 别名]
# 或者: from astropy.table.QTable import add_index [as 别名]
_temperature = [15.513, 15.48, 15.36, 14.404,
                13.37, 12.25, 10.53, 9.30, 8.035,
                7.214, 6.461, 5.531, 4.426, 2.981,
                2.035, 0.884, 0.1818, 0.005770] * u.MK

# density -  g cm^-3
_density = [147.74, 146.66, 142.73, 116.10, 93.35,
            72.73, 48.19, 34.28, 21.958, 15.157,
            10.157, 5.566, 2.259, 0.4483, 0.1528,
            0.042, 0.00361, 1.99e-7] * u.g*u.cm**-3

_d = {'radius': _radius, 'mass': _mass, 'luminosity': _luminosity,
      'temperature': _temperature, 'density': _density}
interior = QTable(_d)
interior.source = 'Turck-Chieze et al. (1988)'
interior.add_index('radius')

# time -  10^9 years
_time = [0, 0.143, 0.856, 1.863, 2.193, 3.020,
         3.977, 4.587, 5.506, 6.074, 6.577, 7.027,
         7.728, 8.258, 8.7566, 9.805] * u.Gyr

# luminosity -  L_sun
_tluminosity = [0.7688, 0.7248, 0.7621, 0.8156,
                0.8352, 0.8855, 0.9522, 1.0, 1.079,
                1.133, 1.186, 1.238, 1.318, 1.399,
                1.494, 1.760] * u.Lsun

# radius -  R_sun
_tradius = [0.872, 0.885, 0.902, 0.924, 0.932,
            0.953, 0.981, 1.0, 1.035, 1.059, 1.082,

# 需要导入模块: from astropy.table import QTable [as 别名]
# 或者: from astropy.table.QTable import add_index [as 别名]

#.........这里部分代码省略.........
                        zip(pt1, pt2)])
        vec /= np.linalg.norm(vec)

        per_vec = perpendicular(vec)

        line_pts = find_path_ends(skel_pts[i], max_pixel, per_vec)

        left_profile, left_dists = \
            profile_line(image, skel_pts[i], line_pts[0])
        right_profile, right_dists = \
            profile_line(image, skel_pts[i], line_pts[1])

        total_profile = np.append(left_profile[::-1], right_profile) * \
            bright_unit

        if noise is not None:
            left_profile, _ = \
                profile_line(noise, skel_pts[i], line_pts[0])
            right_profile, _ = \
                profile_line(noise, skel_pts[i], line_pts[1])
            noise_profile = np.append(left_profile[::-1], right_profile) * \
                bright_unit
        else:
            noise_profile = None

        if distance is not None:
            total_dists = np.append(-left_dists[::-1], right_dists) \
                * u.pix * phys_per_pix
        else:
            total_dists = np.append(-left_dists[::-1], right_dists) \
                * u.pix * deg_per_pix

        if noise is not None:
            if len(total_profile) != len(noise_profile):
                raise ValueError("Intensity and noise profile lengths do not"
                                 " match. Have you applied the same mask to"
                                 " both?")

        line_profiles.append(total_profile)
        line_distances.append(total_dists)
        profile_extents.append([line_pts[0], skel_pts[i], line_pts[1]])

        # Now fit!
        profile_fit, profile_fit_err, red_chisq = \
            gauss_fit(total_dists.value, total_profile.value,
                      sigma=noise_profile)

        profile_fits.append(np.hstack([profile_fit, profile_fit_err]))
        red_chisqs.append(red_chisq)

        if verbose:
            p.subplot(121)
            p.imshow(image, origin='lower')
            p.contour(skeleton, colors='r')
            p.plot(skel_pts[i][1], skel_pts[i][0], 'bD')
            p.plot(line_pts[0][1], line_pts[0][0], 'bD')
            p.plot(line_pts[1][1], line_pts[1][0], 'bD')

            p.subplot(122)
            p.plot(total_dists, total_profile, 'bD')
            pts = np.linspace(total_dists.min().value,
                              total_dists.max().value, 100)
            p.plot(pts, gaussian(pts, *profile_fit), 'r')

            if distance is not None:
                unit = (u.pix * phys_per_pix).unit.to_string()
            else:
                unit = (u.pix * deg_per_pix).unit.to_string()
            p.xlabel("Distance from skeleton (" + unit + ")")
            p.ylabel("Surface Brightness (" + bright_unit.to_string() + ")")
            p.tight_layout()
            p.show()

    profile_fits = np.asarray(profile_fits)
    red_chisqs = np.asarray(red_chisqs)

    # Create an astropy table of the fit results
    param_names = ["Amplitude", "Std Dev", "Background"]
    param_errs = [par + " Error" for par in param_names]
    colnames = param_names + param_errs
    in_bright_units = [True, False, True] * 2
    tab = QTable()

    tab["Number"] = np.arange(profile_fits.shape[0])
    tab.add_index("Number")

    tab["Red Chisq"] = red_chisqs

    for i, (name, is_bright) in enumerate(zip(colnames, in_bright_units)):
        if is_bright:
            col_unit = bright_unit
        else:
            if distance is not None:
                col_unit = (u.pix * phys_per_pix).unit
            else:
                col_unit = (u.pix * deg_per_pix).unit

        tab[name] = profile_fits[:, i] * col_unit

    return line_distances, line_profiles, profile_extents, tab