Merge pull request #23 from Johannes-Sahlmann/plotting-tweaks

Plotting tweaks

pystrometry/gaia_astrometry.py

@@ -90,7 +90,6 @@ def load_data(self, verbose=0, write_modified_iad=0, use_modified_epoch_data=Fal
         if filter_on_frame_uncertainty:
             rse_filter_threshold = 10
             rse = robust_scatter_estimate(self.epoch_data['errda_mas_obs'])
-            # print(rse)
             remove_index = np.where(np.abs(self.epoch_data['errda_mas_obs']) > rse_filter_threshold*rse)[0]
             self.n_filtered_frames += len(remove_index)
             print('RSE cut removed {} frames (of {})'.format(self.n_filtered_frames, self.n_original_frames))
@@ -498,6 +497,8 @@ def plot_individual_orbit(parameter_dict, iad, mapping_dr3id_to_starname=None,
         P_day = parameter_dict['p{}_period_day'.format(planet_number)]
         m2_kg = pjGet_m2(m1_MS * MS_kg, a_m, P_day)
         m2_MJ = m2_kg / MJ_kg
+        m2_MS = m2_kg / MS_kg
+        mass_function = alpha_mas**3 * (u.year.to(u.day))**2 / (P_day**2 * absolute_parallax_mas**3)
 
         # create dictionary that defines the OrbitSystem
         attribute_dict = {
@@ -567,7 +568,9 @@ def plot_individual_orbit(parameter_dict, iad, mapping_dr3id_to_starname=None,
             orbit_descr += '$i={0[p1_incl_deg]:2.{prec}f}\\pm{0[sigma_p1_incl_deg]:2.{prec}f}$ deg\n'.format(dict(parameter_dict), prec=2)
             orbit_descr += '$\\omega={0[p1_omega_deg]:2.{prec}f}\\pm{0[sigma_p1_omega_deg]:2.{prec}f}$ deg\n'.format(dict(parameter_dict), prec=2)
             orbit_descr += '$\\Omega={0[p1_OMEGA_deg]:2.{prec}f}\\pm{0[sigma_p1_OMEGA_deg]:2.{prec}f}$ deg\n'.format(dict(parameter_dict), prec=2)
+            orbit_descr += '$f(m)={0:2.{prec}f}$ Msun\n'.format(mass_function, prec=2)
             orbit_descr += '$M_1={0:2.{prec}f}$ Msun\n'.format(m1_MS, prec=2)
+            orbit_descr += '$M_2={0:2.{prec}f}$ Msun\n'.format(m2_MS, prec=2)
             orbit_descr += '$M_2={0:2.{prec}f}$ Mjup\n'.format(m2_MJ, prec=2)
 
         else:
@@ -620,6 +623,10 @@ def plot_individual_orbit(parameter_dict, iad, mapping_dr3id_to_starname=None,
         axp.title = '{} ({})'.format(iad.source_id, mag_str)
         name_seed = '{}'.format(iad.source_id)
 
+    if ('title' in parameter_dict.keys()) and (parameter_dict['title'] is not None):
+        axp.title = parameter_dict['title']
+
+
     argument_dict = {'plot_dir'             : plot_dir, 'ppm_panel': True,
                      'frame_residual_panel' : True, 'orbit_only_panel': True,
                      'ppm_description'      : ppm_description,
@@ -651,6 +658,7 @@ def plot_individual_orbit(parameter_dict, iad, mapping_dr3id_to_starname=None,
 # small helper dictionary to become tolerant against different variable naming, e.g.
 # internal DPAC naming versus Gaia archive naming
 astrometric_parameter_name_equivalences = {'varpi': 'parallax',
+                                           'varpi': 'plx_mas',
                                            'plx_mas': 'parallax',
                                            'plx_mas': 'absolute_plx_mas',
                                            'absolute_plx_mas': 'parallax',
@@ -784,7 +792,10 @@ def plot_individual_ppm(parameter_dict, iad, plot_dir=os.path.expanduser('~')):
 
     source_id = np.int64(parameter_dict['sourceId'])
     # axp.title = 'Gaia DR4 {} ({}, {})'.format(source_id, 'LPC', mag_str)
-    axp.title = 'Gaia {} ({})'.format(source_id, mag_str)
+    if ('title' in parameter_dict.keys()) and (parameter_dict['title'] is not None):
+        axp.title = parameter_dict['title']
+    else:
+        axp.title = 'Gaia {} ({})'.format(source_id, mag_str)
 
     name_seed = parameter_dict.get('name_seed', 'Gaia_{}'.format(source_id))
 

pystrometry/pystrometry.py

@@ -2225,7 +2225,7 @@ def __init__(self, attribute_dict=None):
 
         if len(outlier_1D_index) != 0:
             print('MJD of outliers:')
-            for ii in np.unique(outlier_1D_index.astype(np.int)):
+            for ii in np.unique(outlier_1D_index.astype(int)):
                 print('{:.12f}'.format(T['MJD'][ii]), end=',')
             print()
 
@@ -2468,7 +2468,7 @@ def plot(self, argument_dict=None):
                                     '$\\bar\\sigma_\Lambda$={:2.3f} mas'.format(
                 len(self.data.epoch_data), self.data.n_original_frames, np.std(self.residuals), np.mean(self.data.epoch_data['sigma_da_mas']))
             if 'excess_noise' in argument_dict.keys():
-                argument_dict['frame_omc_description'] += '\nexN = {:2.2f}, mF = {:2.0f}'.format(
+                argument_dict['frame_omc_description'] += '\nexN = {:2.2f} mas, mF = {:2.0f}'.format(
             argument_dict['excess_noise'], argument_dict['merit_function'])
 
         if argument_dict['orbit_signal_description'] == 'default':
@@ -3380,7 +3380,7 @@ def set_binned_quantities(self):
 
         if len(outlier_1D_index) != 0:
             print('MJD of outliers:')
-            for ii in np.unique(outlier_1D_index.astype(np.int)):
+            for ii in np.unique(outlier_1D_index.astype(int)):
                 print('{:.12f}'.format(T['MJD'][ii]), end=',')
             print()
 

pystrometry/utils/archives.py

@@ -196,10 +196,14 @@ def query_dpcg_epochastrometry(connection, out_dir, tag='dpcgdata', query=None,
             from mdb_gaia_starobject_088 so
             join dr3_ops_cs36_mv.dgdreq58_rejected_cu4transitids_astro transrej using (sourceid)
             join lateral unnest(filter_transits(transits,rejected_cu4transitids_astro)) as t on true
-            -- provide source id list (comma separated)
-            where sourceid in ({selected_source_id_string})
             """
 
+            if selected_source_id_string != None:
+                query += f"""-- provide source id list (comma separated) 
+                where sourceid in ({selected_source_id_string})"""
+            else:
+                query += """ where sourceid != 0"""
+
         dpcg_df = pd.read_sql(query, connection)
         dpcg_df.to_parquet(out_file)
         logging.info(f'Wrote {len(dpcg_df)} rows to {out_file}')

pystrometry/utils/du437_tools.py

@@ -162,7 +162,7 @@ def draw_monte_carlo_solution_parameters(row, n_mc=1000, parameters=None, indice
 
 
         cov_matrix = covariance_matrix(row, parameters=parameters)
-        solution_parameters = row[parameters].astype(np.float)
+        solution_parameters = row[parameters].astype(float)
 
         # draw random data accounting for covariances
         if np.isnan(cov_matrix).any(): # this corresponds to circular orbits
@@ -327,7 +327,7 @@ def sample_solution_parameters_monte_carlo(self, index, n_mc=1000):
         #         # print(cov_matrix)
         #     else:
         cov_matrix = covariance_matrix(self._obj.loc[index], parameters=parameters)
-        solution_parameters = self._obj.loc[index, parameters].astype(np.float)
+        solution_parameters = self._obj.loc[index, parameters].astype(float)
 
         # draw random data accounting for covariances
         if np.isnan(cov_matrix).any(): # this corresponds to circular orbits
@@ -400,7 +400,7 @@ def compute_ge_mc(solution, n_mc=1000):
                 # print(cov_matrix)
             else:
                 cov_matrix = covariance_matrix(solution)[5:9, 5:9]
-            thiele_innes_parameters = solution[parameters].astype(np.float)
+            thiele_innes_parameters = solution[parameters].astype(float)
 
             # draw random data accounting for covariances
             if np.isnan(cov_matrix).any(): # this corresponds to circular orbits
@@ -543,7 +543,7 @@ def plot_cmd(self, label=None, title=None, ax=None, show_cutoff=False, **kwargs)
             ax = pl.gca()
 
         self._obj.plot('bp_rp', 'absolute_phot_g_mean_mag', kind='scatter', ax=ax, label=label, **kwargs)  # , c='parallax'
-        pl.title('{} ({} sources)'.format(title, len(self._obj)))
+        ax.set_title('{} ({} sources)'.format(title, len(self._obj)))
         ax.invert_yaxis()
         if show_cutoff:
             colour_cutoff = 3.5
@@ -557,10 +557,10 @@ def cutoff(x):
             pl.plot(x, cutoff(x), 'k-')
             pl.plot([colour_cutoff, colour_cutoff], [cutoff(colour_cutoff), 0], 'k-')
 
-        pl.xlabel('$G_\mathrm{BP}-G_\mathrm{RP}$')
-        pl.ylabel('$M_G$')
+        ax.set_xlabel('$G_\mathrm{BP}-G_\mathrm{RP}$')
+        ax.set_ylabel('$M_G$')
 
-        pl.show()
+        # pl.show()
         # return ax
 
 @pd.api.extensions.register_dataframe_accessor("nssmc")

setup.cfg

@@ -40,10 +40,10 @@ github_project = https://github.com/Johannes-Sahlmann/pystrometry
 # install_requires should be formatted as a comma-separated list, e.g.:
 # install_requires = astropy, scipy, matplotlib
 
-install_requires = astropy, linearfit, matplotlib, scipy, astroquery, sympy
+install_requires = astropy, linearfit>=1.0.2, matplotlib, scipy, astroquery, sympy
 
 # version should be PEP440 compatible (https://www.python.org/dev/peps/pep-0440/)
-version = 0.4.3
+version = 0.5.0
 # Note: you will also need to change this in your package's __init__.py
 minimum_python_version = 3.7