dfm / tess-atlas

MIT License
9 stars 8 forks source link

Initial lc phase plot error #228

Open avivajpeyi opened 2 years ago

avivajpeyi commented 2 years ago
cat july12_cat/log_pe/pe_28112380_778.log
WARNING (theano.tensor.blas): Using NumPy C-API based implementation for BLAS functions.
[00:00:05 - TESS-ATLAS-RUNNER] run_toi(5581)
[00:00:05 - TESS-ATLAS-RUNNER] Executing july12_cat/0.2.1.dev64+gc7fa3a0/toi_5581.ipynb
[00:03:42 - TESS-ATLAS-RUNNER] Preprocessing july12_cat/0.2.1.dev64+gc7fa3a0/toi_5581.ipynb failed:

 An error occurred while executing the following cell:
------------------
plot_phase(
    tic_entry,
    planet_transit_model,
    initial_params=init_params,
    plot_all_datapoints=True,
)
------------------

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
/tmp/ipykernel_8468/1078644873.py in <module>
----> 1 plot_phase(
      2     tic_entry,
      3     planet_transit_model,
      4     initial_params=init_params,
      5     plot_all_datapoints=True,

/fred/oz200/avajpeyi/projects/tess-atlas/src/tess_atlas/plotting/matplotlib_plots.py in plot_phase(tic_entry, model, inference_data, initial_params, data_bins, plot_error_bars, plot_data_ci, plot_all_datapoints, zoom_y_axis, plot_label, num_lc)
    231         # Plot the folded lightcurve model
    232         if len(initial_lightcurves) > 0:
--> 233             init_x, init_y = fold_lightcurve_models(
    234                 lc,
    235                 initial_lightcurves[..., i],

/fred/oz200/avajpeyi/projects/tess-atlas/src/tess_atlas/plotting/plotting_utils.py in fold_lightcurve_models(lc, lc_models, t0s, periods, plt_min, plt_max)
     46     # interpolate each lightcurve model to have same lengths
     47     for i, (x, y) in enumerate(zip(xs, ys)):
---> 48         f = interp1d(x, y, fill_value="extrapolate")
     49         ys[i] = f(longest_x)
     50     return longest_x, np.array(ys)

/fred/oz200/avajpeyi/envs/tess/lib/python3.8/site-packages/scipy/interpolate/interpolate.py in __init__(self, x, y, kind, axis, copy, bounds_error, fill_value, assume_sorted)
    472         # Interpolation goes internally along the first axis
    473         self.y = y
--> 474         self._y = self._reshape_yi(self.y)
    475         self.x = x
    476         del y, x  # clean up namespace to prevent misuse; use attributes

/fred/oz200/avajpeyi/envs/tess/lib/python3.8/site-packages/scipy/interpolate/polyint.py in _reshape_yi(self, yi, check)
    108                                            self._y_extra_shape[:-self._y_axis])
    109             raise ValueError("Data must be of shape %s" % ok_shape)
--> 110         return yi.reshape((yi.shape[0], -1))
    111
    112     def _set_yi(self, yi, xi=None, axis=None):

ValueError: cannot reshape array of size 0 into shape (0,newaxis)
ValueError: cannot reshape array of size 0 into shape (0,newaxis)

[00:03:42 - TESS-ATLAS-RUNNER] TOI 5581 execution complete: False (216.73s)

Other TOIs with the same issue:

359
710
978
1192
1894
1895
2005
2151
2168
2299
2321
2326
2338
2485
2537
2670
3487
4127
4144
4344
4348
4355
4356
4358
4555
4562
4635
4706
5087
5149
5153
5235
5396
5542
5571
5581
5619
5624
avivajpeyi commented 2 years ago
import pandas as pd
x = "359 710 978 1192 1894 1895 2005 2151 2168 2299 2321 2326 2338 2485 2537 2670 3487 4127 4144 4344 4348 4355 4356 4358 4555 4562 4635 4706 5087 5149 5153 5235 5396 5542 5571 5581 5619 5624".split()
p = [f"toi_{i}_files/tic_data.csv" for i in x]
df = pd.concat([pd.read_csv(i) for i in p])
df[['TOI int', 'Single Transit', 'Multiplanet System', 'Planet SNR']]
   TOI int  Single Transit  Multiplanet System  Planet SNR
0      359           False               False   10.000000
0      710           False               False   33.773872
0      978           False               False   32.000000
0     1192           False               False   42.492947
0     1894           False               False   24.061493
0     1895           False               False   40.682007
0     2005           False               False   36.000000
0     2151           False               False   43.000000
0     2168           False               False   68.000000
0     2299           False               False    7.914845
0     2321            True               False    7.323930
0     2326           False               False   21.000000
0     2338           False               False   55.000000
0     2485           False               False   49.000000
0     2537           False               False   70.000000
0     2670           False               False   21.000000
0     3487           False               False   37.000000
0     4127           False               False   47.000000
0     4144           False               False   87.000000
0     4344           False               False    7.574218
0     4348           False               False    7.539507
0     4355           False               False    7.247839
0     4356           False               False    7.153257
0     4358           False               False   11.755434
0     4555           False               False    7.200012
0     4562           False               False   31.857893
0     4635           False               False   90.265900
0     4706           False               False   50.000000
0     5087           False               False   23.000000
0     5149           False               False  158.000000
0     5153           False               False   34.000000
0     5235           False               False   11.000000
0     5396           False               False   14.000000
0     5542           False               False   23.000000
0     5571           False               False   31.000000
0     5581           False               False   43.000000
0     5619           False               False   10.000000
0     5624           False               False   19.000000
avivajpeyi commented 2 years ago

Ok! all of these appear to have "flat-lined" inital lightcurves --> period probs a nan

Screen Shot 2022-07-15 at 12 23 15 pm Screen Shot 2022-07-15 at 12 23 20 pm

Yep! tmin == tmax!

Screen Shot 2022-07-15 at 12 25 42 pm
avivajpeyi commented 2 years ago

It seems like there are some single-transit systems which have non-zero period estimates on ExoFop.

We set the system to be a single transit system if the following holds true: https://github.com/dfm/tess-atlas/blob/8427171e9be11dc8f5839f9aff99de3b5294b758/src/tess_atlas/data/planet_candidate.py#L49

If the system is not a single-transit system, then our model uses tmin, tmax, where:

tmin: https://github.com/dfm/tess-atlas/blob/8427171e9be11dc8f5839f9aff99de3b5294b758/src/tess_atlas/data/planet_candidate.py#L90

tmax: https://github.com/dfm/tess-atlas/blob/8427171e9be11dc8f5839f9aff99de3b5294b758/src/tess_atlas/data/planet_candidate.py#L83

If not, we use

lc_5542

lc_5571

lc_5581

lc_5619

lc_5624

lc_359

lc_710

lc_978

lc_1192

lc_1894

lc_1895

lc_2005

lc_2151

lc_2168

lc_2299

lc_2321

lc_2326

lc_2338

lc_2485

lc_2537

lc_2670

lc_3487

lc_4127

lc_4144

lc_4344

lc_4348

lc_4355

lc_4356

lc_4358

lc_4555

lc_4562

lc_4635

lc_4706

lc_5087

lc_5149

lc_5153

lc_5235

lc_5396

import matplotlib.pyplot as plt
from tess_atlas.data import TICEntry
from tess_atlas.plotting import (
    plot_raw_lightcurve
)

from tqdm.auto import tqdm
import os

FAILED_TOIS = "359 710 978 1192 1894 1895 2005 2151 2168 2299 2321 2326 2338 2485 2537 2670 3487 4127 4144 4344 4348 4355 4356 4358 4555 4562 4635 4706 5087 5149 5153 5235 5396 5542 5571 5581 5619 5624".split()
FAILED_TOIS = [int(i) for i in FAILED_TOIS]

OUT = "out"
os.makedirs(OUT, exist_ok=True)

for toi_n in tqdm(FAILED_TOIS, desc="Plotting TOI"):
    tic = TICEntry.load(toi_n)
    fig = plot_raw_lightcurve(tic, save=False)
    fig.savefig(f"{OUT}/lc_{toi_n}.png")
    plt.close(fig)
avivajpeyi commented 2 years ago

Taking a closer look at TOI 359, i think that some of the t0 may be off? Also... how was T=16 estimated for this system lol

Here i have plotted [t0, t0+T, t0+2T, .... ]. Seems like the exofop t0 isnt a valid epoch?

test

hmm idk whats going on with this system.. bls

avivajpeyi commented 2 years ago

Maybe situations when tmin==tmax means that the planet only has a single transit?

This works for TOI 5153

Screen Shot 2022-07-15 at 5 13 27 pm

Here, tmin==tmax --> and this looks like a single-transit system

However, this does not work for TOI 4635

Screen Shot 2022-07-15 at 5 13 40 pm

Here, tmin==tmax , but there clearly is a second eclipse @ t0 + T

avivajpeyi commented 2 years ago

@dfm -- any suggestions on this issue?

dfm commented 2 years ago

Good question. I think this'll have something to do with floor vs. ceil and ±1s in our definitions of tmin and tmax. It should be possible to define them such that tmin == tmax only when you have a single transit. Can you try experimenting with that to see if you can come up with a consistent definition that works in these test cases (one that's a single transit and one that's a double)?

avivajpeyi commented 2 years ago

@dfm -- im having trouble getting reasonable epoch values for some of these failed TOIs.

I think the ExoFop epoch may be slightly off?

Eg see the below

Screen Shot 2022-07-19 at 4 23 52 pm
%pip install lightkurve -q
import lightkurve as lk
import pandas as pd

EXOFOP_TOI = "https://exofop.ipac.caltech.edu/tess/download_toi.php?output=csv&toi={}"  

BJD = 2457000
TOI = 1895

class ToiData:

    def __init__(self, toi):
        self.toi = toi
        self.lc = lk.search_lightcurve(f"TOI {toi}", author="SPOC").download_all().stitch().remove_nans() -1
        self.params = pd.read_csv(EXOFOP_TOI.format(toi)).to_dict('records')[0]
        self.T = self.params['Period (days)']
        t0 = self.params['Epoch (BJD)']-BJD
        self.epochs = [t0 - self.T, t0, t0 + self.T]

    def plot(self):
        ax = self.lc.scatter()
        ax.plot(self.epochs, [0,0,0], 'ro-')
        ax.vlines([self.epochs[1]], min(self.lc.flux), 0, 'r', label='Exofop Epoch')
        ax.vlines(self.epochs, min(self.lc.flux), 0, 'r', ls='--', label='Epoch+/-period')
        ax.legend()

ToiData(TOI).plot()
avivajpeyi commented 2 years ago

Ok maybe i shouldn't test with these TOIs where the epoch might be off (eg 1895, 359). Maybe these are just some TOIs that we wont be able to analyse.

dfm commented 2 years ago

If the TOI list has incorrect entries then I think it's totally fine to punt on that for now!