#!/usr/bin/perl

=head1 NAME

lookup_mpc - find detections from MPC obs format in current MOPS DB

=head1 SYNOPSIS

lookup_mpc < OBSFILE

  OBSFILE : MPC-formatted observations
  --debug : dump data, don't insert anything
  --help : show this manpage

=head1 DESCRIPTION

=cut


use strict;
use warnings;
use Data::Dumper;
use Getopt::Long;
use Pod::Usage;
use Mail::Send;
use Astro::Time;

# Just a stub for quick debugging.
use PS::MOPS::Lib qw(:all);
use PS::MOPS::MPC;
use PS::MOPS::DC::Instance;
use PS::MOPS::DC::Field;
use PS::MOPS::DC::Tracklet;
use PS::MOPS::DC::DerivedObject;


# Detection statuses.
our $SUBMISSION_PENDING_NEO = 'P';
our $SUBMISSION_PENDING_STD = 'Q';
our $SUBMISSION_NEO = 'S';
our $SUBMISSION_STD = 'T';
our $SUBMISSION_OUTREACH = 'O';

# Start program here.
my $instance_name;
my $debug;
my $help;
GetOptions(
    'instance=s' => \$instance_name,
    debug => \$debug,
    help => \$help,
) or pod2usage(2);
pod2usage(-verbose => 3) if $help;


my $inst = PS::MOPS::DC::Instance->new(DBNAME => $instance_name);
my $dbname = $inst->dbname;
my $logger = $inst->getLogger();
my $dbh = $inst->dbh;
my $lookup_sth = $dbh->prepare(<<"SQL") or die $dbh->errstr;
/* args are P1.dec, P2.dec, P1.dec, P2.dec, P1.ra, P2.ra */
select det_id 
from detections d join fields f using(field_id)
where f.fpa_id=?
and abs(
    degrees(
        /* spherical distance from tracklet to target field center */
        acos(least(1.0, /* handle slight roundoff error causing arg > 1.0 */
            sin(radians(d.dec_deg)) * sin(radians(?))
            + cos(radians(d.dec_deg)) * cos(radians(?)) * cos(radians(d.ra_deg - ?))
        ))
    )
) < 0.000555 /* 2 arcsec */
SQL


my %already_inserted;
my %trk_to_insert;


# First read all the MPC-formatted detections and arrange into 
# designations (tracklets).
my $tracklets;
$tracklets = read_mpc();        # read all from STDIN
while (my ($desig, $detref) = each %{$tracklets}) {
    # 1. Check that none of the tracklet's detections have been submitted.
    # 2. Set up insert data.
    # 3. Insert into table & retrieve designation.
    # 4. Submit via email.

    # This is a tracklet; fetch it and add the detections.
    my $can_insert = 1;
    foreach my $det (@{$detref}) {
        # Look up det in DB.
        my $found_desig;
        my $det_id;
        if ($det_id = lookup_det($det)) {
            emit(sprintf "Found %s as %s/%s.\n", $desig, $dbname, $det_id);
        }
        else {
            emit(sprintf "$desig not found at RA=%.6f DEC=%.6f in $dbname/%s.\n", $det->{ra}, $det->{dec}, $det->{fpa_id});
        }
    }
}

if ($debug) {
    # Sanity check.
    foreach my $desig (keys %trk_to_insert) {
        print Dumper(\%trk_to_insert);
    }
}
exit;


sub emit {
    print @_;
}


sub read_mpc {
    my $line;
    my %desigs;
    while ($line = <STDIN>) {
        chomp $line;
        next if !$line or $line =~ /^\s*$/ or $line =~ /^[A-Z]/;
        my $foo = PS::MOPS::MPC::mpc_split($line);
        my $desig = $foo->{DESIGNATION} or die "didn't get designation from $line";
        die "Not a PS1 observation: $foo->{OBSERVATORY}\n" unless $foo->{OBSERVATORY} eq 'F51';

        my $f = lookup_field($foo->{DATE})
            or die "can't find matching field for $foo->{DATE} in $dbname\n";
        my $det = {
            desig => $foo->{DESIGNATION},
            epoch_mjd => $f->{epoch_mjd},
            fpa_id => $f->{fpa_id},
            survey_mode => $f->{survey_mode},
            ra => str2deg($foo->{RA}, 'H'),
            dec => str2deg($foo->{DEC}, 'D'),
            mag => $foo->{MAG},
            filter_id => $foo->{BAND},
            obscode => $foo->{OBSERVATORY},
        };
        push @{$desigs{$desig}}, $det;
    }

    return \%desigs;
}


sub lookup_det {
    my ($det) = @_;
    $lookup_sth->execute(
        $det->{fpa_id}, $det->{dec}, $det->{dec}, $det->{ra}
    ) or die $lookup_sth->errstr;
    my ($near) = $lookup_sth->fetchrow_array();
    return $near;
}


sub lookup_field {
    # Given an MPC date string, return an FPA id for the nearest PS1 
    # exposure within +/- 1 sec.  Date format is
    # 2011 01 23.23452
    my ($date) = @_;
    $date =~ s/^\s+|\s+$//g;    # strip
    my ($y, $m, $d) = split /\s+/, $date;
    my $df = ($d - int($d));
    $d = int($d);
    my $mjd = cal2mjd($d, $m, $y, $df);
    my $field_i = modcf_retrieve($inst, 
        begin_mjd => $mjd - 5/100_000,
        end_mjd => $mjd + 5/100_000,
    );

    my $field = $field_i->next;
    if ($field) {
        return {
            epoch_mjd => $field->epoch,
            fpa_id => $field->fpaId,
            survey_mode => $field->surveyMode,
        };
    }
    else {
        return undef;
    }
}


sub read_tracklet {
    my ($desig, $tid) = @_;
    my $tdata = {};
    my $trk = modct_retrieve($inst, trackletId => $tid) or die "can't retrieve tracklet ID $tid\n";
    foreach my $det (sort {$a->epoch <=> $b->epoch} @{$trk->detections}) {
        my $field = modcf_retrieve($inst, fieldId => $det->fieldId) or die "can't retrieve field ID $det->fieldId\n"; 

        push @{$tdata->{$desig}}, {
            desig => $desig,
            det_id => $det->detId,
            epoch_mjd => $det->epoch,
            fpa_id => $field->fpaId,
            survey_mode => $field->{surveyMode},
            ra => $det->ra,
            dec => $det->dec,
            mag => $det->mag,
            filter_id => $det->filter,
            obscode => $det->obscode,
        };
    }

    return $tdata;
}
