Subs::Attitude (version $)


package Subs::Attitude;
##############################################################################
#
# DESCRIPTION: Extract attitude data. Swift attitude data appears
# DESCRIPTION: in several places in the telemetry. This subroutine
# DESCRIPTION: collects all the attitude data into a single quaternion-
# DESCRIPTION: based attitude file and removed redundant information.
# DESCRIPTION:
# DESCRIPTION: The main source of attitude data is the "ACS packets"
# DESCRIPTION: (APID 404). The BAT also places ACS records in the
# DESCRIPTION: headers of each LDP and in the body of a special
# DESCRIPTION: LDP generated suring slews.
# DESCRIPTION: The XRT places ACS records in the header for each CCD
# DESCRIPTION: frame, however these are represented by "floats" instead
# DESCRIPTION: of "doubles", so they are not used if higher accuracy
# DESCRIPTION: data are available.
# DESCRIPTION:
# DESCRIPTION: The sub checkDatabaseForNominalPointing reads the database
# DESCRIPTION: file "/sw/nominal_pnt.rdb" from the directory indicated by
# DESCRIPTION: env var $RDB_TABLES_DIR (normally set in aps.config;
# DESCRIPTION: defaults to /aps/db/rdb/tables).  If the sequence is found
# DESCRIPTION: in this table, those values are used for the nominal
# DESCRIPTION: pointing; this allows for manual override when the pointing
# DESCRIPTION: cannot be reliably determined automatically.  Otherwise, the
# DESCRIPTION: mean pointing is calculated and compared with the requested
# DESCRIPTION: pointing (SACSREQ, read from the engineering housekeeping
# DESCRIPTION: file), if these disagree an E2 error is logged and an email
# DESCRIPTION: notification is sent.  Currently the tolerance is 10 arcmin.
# DESCRIPTION: Someone should then review the sequence and add the correct
# DESCRIPTION: position to nominal_pnt.rdb, then rerun the sequence.
# DESCRIPTION:
# DESCRIPTION: Here is a pseudocode outline of the algorithm used to
# DESCRIPTION: determine the nominal attitude:
# DESCRIPTION:
# DESCRIPTION:  IF obs in ADB [nominal pointing database]
# DESCRIPTION:     pnt_db = ADB[obs]
# DESCRIPTION:  ELSE
# DESCRIPTION:     pnt_sacs = determineSACSREQ(sw<obs>sen.hk)
# DESCRIPTION:     # 1-2% of the time, determineSACSREQ fails...
# DESCRIPTION:     pnt_aspect = aspect(sw<obs>[pu]at.fits)
# DESCRIPTION:     IF pnt_sacs AND pnt_aspect THEN
# DESCRIPTION:        IF separation(pnt_sacs, pnt_aspect) > 10 arcmin THEN
# DESCRIPTION:           E2 ERROR SACSREQ and aspect positions inconsistent
# DESCRIPTION:           SEND EMAIL NOTIFICATION
# DESCRIPTION:        ENDIF
# DESCRIPTION:        IF minRollDelta(pnt_sacs, pnt_aspect) > 2 deg THEN
# DESCRIPTION:           E1 ERROR
# DESCRIPTION:        ENDIF
# DESCRIPTION:     ELSEIF pnt_aspect
# DESCRIPTION:        IF pnt_aspect has ERROR
# DESCRIPTION:           E2 ERROR Unable to determine nominal pointing
# DESCRIPTION:           SEND EMAIL NOTIFICATION
# DESCRIPTION:        ENDIF
# DESCRIPTION:     ENDIF
# DESCRIPTION:  ENDIF
# DESCRIPTION:
#
# HISTORY: $Log: Attitude.pm,v $
# HISTORY: Revision 1.63  2016/11/04 13:04:49  apsop
# HISTORY: Add UTCFINIT to attitude files so prefilter can use it to address clock drift. Pass prefilter jump corrected (.pat) instead of uncorrected (.sat) attitude if available.
# HISTORY:
# HISTORY: Revision 1.62  2014/08/14 09:19:47  apsop
# HISTORY: Send an email notification to watchers if nominal pointing cannot be
# HISTORY: determined by any method.  Include job title in both nominal-pointing
# HISTORY: emails.
# HISTORY:
# HISTORY: Revision 1.61  2014/02/27 08:26:31  apsop
# HISTORY: ATTFLAG comment explains value.
# HISTORY: Sequence and seqprocnum included on subject line of
# HISTORY: aspect/SACSREQ Inconsistent email.
# HISTORY:
# HISTORY: Revision 1.60  2013/04/29 22:08:50  apsop
# HISTORY: Add seqprocnum to the email sent out if aspect & sacsreq are inconsistent
# HISTORY:
# HISTORY: Revision 1.59  2013/03/28 09:00:58  apsop
# HISTORY: Added more comments.
# HISTORY:
# HISTORY: Revision 1.58  2013/03/14 10:02:11  apsop
# HISTORY: Math.pm has been added to Util, so use Util::Math
# HISTORY:
# HISTORY: Revision 1.57  2012/09/08 07:03:56  apsop
# HISTORY: Added object, sequence, target to email alert for inconsitent
# HISTORY: SACSREQ/aspect positions.
# HISTORY:
# HISTORY: Revision 1.56  2012/09/05 03:44:12  apsop
# HISTORY: Modified command that reads nominal_pnt.rdb to not return header.
# HISTORY:
# HISTORY: Revision 1.55  2012/08/15 06:53:37  apsop
# HISTORY: Check database nominal_pnt.rdb for manual pointing override;
# HISTORY: compare mean aspect solution with requested pointing (SACSREQ);
# HISTORY: add time keywords (TSTART, TSTOP, DATE, etc.) to all extensions.
# HISTORY: [Bob W.]
# HISTORY:
# HISTORY: Revision 1.54  2007/10/10 20:32:29  apsop
# HISTORY: Put tolerance for target-pointing GTI back to 1.0 .
# HISTORY:
# HISTORY: Revision 1.53  2007/10/08 21:03:56  apsop
# HISTORY: For targ pointing GTI, reduced tolerance to 0.1 degrees, and rename STDGTI extension to GTI.
# HISTORY:
# HISTORY: Revision 1.52  2007/09/28 17:14:48  apsop
# HISTORY: Change writing of ATTFLAG keyword so that it is a string rather than an integer.
# HISTORY:
# HISTORY: Revision 1.51  2007/09/17 15:32:15  apsop
# HISTORY: Set chatter for attjumpcorr to 3.
# HISTORY:
# HISTORY: Revision 1.50  2007/09/11 18:22:23  apsop
# HISTORY: Turn up chatter on attjumpcorr tasked to 5.
# HISTORY:
# HISTORY: Revision 1.49  2007/04/01 19:08:48  apsop
# HISTORY: Set chatter for attjumpcorr to 3.
# HISTORY:
# HISTORY: Revision 1.48  2007/02/22 18:35:33  apsop
# HISTORY: Fix subtle, intemittent bug in running acs2fits w/ setting $_ .
# HISTORY:
# HISTORY: Revision 1.47  2007/02/08 21:33:12  apsop
# HISTORY: Write ATTFLAG keyword to all extensions.
# HISTORY:
# HISTORY: Revision 1.46  2007/01/31 20:22:14  apsop
# HISTORY: Do attitude correction with attjumpcorr and set attitude correction flags.
# HISTORY:
# HISTORY: Revision 1.45  2006/10/01 18:55:57  apsop
# HISTORY: Insert date keywords into the attitude file so that new version of aspect task will work properly.
# HISTORY:
# HISTORY: Revision 1.44  2006/06/28 18:54:42  apsop
# HISTORY: Create new target pointing GTI file.
# HISTORY:
# HISTORY: Revision 1.43  2006/06/15 22:11:52  apsop
# HISTORY: Make pointing and slew GTIs more restrictive, so that times beyond end measurements are excluded.
# HISTORY:
# HISTORY: Revision 1.42  2006/04/27 16:28:12  apsop
# HISTORY: Make a lack of any attitude file a fatal error.
# HISTORY:
# HISTORY: Revision 1.41  2006/03/13 17:39:17  apsop
# HISTORY: Use XRT attitude info in 991 seqs. Fix bug in choosing GTI for mean pointing determination.
# HISTORY:
# HISTORY: Revision 1.40  2006/01/29 19:37:26  apsop
# HISTORY: Add 60 secs of buffer at each end of the range for the att/orbit file.
# HISTORY:
# HISTORY: Revision 1.39  2005/11/20 21:26:10  apsop
# HISTORY: Unexpected OO behaviour for $inter_file requires seperate call to append_to.
# HISTORY:
# HISTORY: Revision 1.38  2005/11/20 20:30:03  apsop
# HISTORY: Check for presence of BUS_V column before making and using ACS_DATA extension in attitude file.
# HISTORY:
# HISTORY: Revision 1.37  2005/11/08 16:58:18  apsop
# HISTORY: Change fdiff to ftdiff and set reltol to 10e-8.  Use caldb to get alignment file in prefilter and abberator.
# HISTORY:
# HISTORY: Revision 1.36  2005/09/26 21:32:35  apsop
# HISTORY: Make not_pointing GTI file.
# HISTORY:
# HISTORY: Revision 1.35  2005/07/29 14:04:03  apsop
# HISTORY: look for ACS packets that are labeled as head2
# HISTORY:
# HISTORY: Revision 1.34  2005/06/01 17:38:51  apsop
# HISTORY: More robust algorithm for selecting GTI file to use for determining the mean pointing.
# HISTORY:
# HISTORY: Revision 1.33  2005/04/19 16:05:48  apsop
# HISTORY: Check for existence of gti files before setting keywords.
# HISTORY:
# HISTORY: Revision 1.32  2005/03/25 21:25:34  apsop
# HISTORY: Fix bug in setting TSTOP for GTIs.  Set gti extname to GTI.
# HISTORY:
# HISTORY: Revision 1.31  2005/03/25 20:18:04  apsop
# HISTORY: Set TSTART and TSTOP for GTI files.
# HISTORY:
# HISTORY: Revision 1.30  2005/03/15 20:04:42  apsop
# HISTORY: Fix bug in call for getting xrt eng hk file name.
# HISTORY:
# HISTORY: Revision 1.30  2005/03/15 19:03:26  apsop
# HISTORY: Check for presence of xrt no position message when making tdrss cat file.
# HISTORY:
# HISTORY: Revision 1.29  2005/02/18 01:53:09  apsop
# HISTORY: Only use secondary attitude data if primary data is not available. Fix up xrt modes.
# HISTORY:
# HISTORY: Revision 1.28  2005/02/14 19:23:39  apsop
# HISTORY: Remove batatt call. Instead pickup bat attitude info from bat2fits output.
# HISTORY:
# HISTORY: Revision 1.27  2005/02/10 02:52:14  apsop
# HISTORY: Get start and stop times from attitude file, as StartEndTimes has not been run yet.
# HISTORY:
# HISTORY: Revision 1.26  2005/02/08 18:24:17  apsop
# HISTORY: Abberation correction using abberator tool with time test.  Move calculation of attorb file to this class as it is needed for abberator.
# HISTORY:
# HISTORY: Revision 1.25  2005/01/12 17:19:13  apsop
# HISTORY: Exclude safeholds from settling and pointing GTIs.  Include 10arcmin flag in pointing GTI.
# HISTORY:
# HISTORY: Revision 1.24  2004/12/10 02:16:09  apsop
# HISTORY: Do a diff of attitude files after normalization, so we have a record of what changed.
# HISTORY:
# HISTORY: Revision 1.23  2004/12/02 18:55:50  apsop
# HISTORY: Normalize quaternions before passing to aspect.
# HISTORY:
# HISTORY: Revision 1.22  2004/11/02 21:12:15  apsop
# HISTORY: Set DATE keyword in attitude file.
# HISTORY:
# HISTORY: Revision 1.21  2004/09/03 00:26:01  apsop
# HISTORY: Temporarily remove bat lc attitude information from attitude file.
# HISTORY:
# HISTORY: Revision 1.20  2004/08/16 15:23:09  apsop
# HISTORY: Turn on writing of HISTORY keywords.
# HISTORY:
# HISTORY: Revision 1.19  2004/07/06 19:59:32  apsop
# HISTORY: Add test for existence of nfi GTI file.
# HISTORY:
# HISTORY: Revision 1.18  2004/06/29 14:34:02  apsop
# HISTORY: New ontarget gti is pointing ANDed with nfi.  Use as input to aspect tool.
# HISTORY:
# HISTORY: Revision 1.17  2004/06/14 14:26:23  apsop
# HISTORY: Write in name of alignment file into attitude file.
# HISTORY:
# HISTORY: Revision 1.16  2004/06/10 19:55:03  apsop
# HISTORY: Add updating of TSTART/TSTOP keywords in second extension.
# HISTORY:
# HISTORY: Revision 1.15  2004/05/28 19:16:07  apsop
# HISTORY: Changes to attitude file format
# HISTORY:
# HISTORY: Revision 1.14  2004/05/06 20:02:33  dah
# HISTORY: Add version number back into the header comments.
# HISTORY:
# HISTORY: Revision 1.13  2004/05/04 16:30:51  dah
# HISTORY: Set ra/dec to 0/0 if not attitude data.  Issue error.
# HISTORY:
# HISTORY: Revision 1.12  2004/04/28 13:45:50  dah
# HISTORY: Change xrt header and engineering file modes to fit new standard.
# HISTORY:
# HISTORY: Revision 1.11  2004/04/16 20:21:18  dah
# HISTORY: Begin using embedded history records
# HISTORY:
#
# VERSION: $Revision: 1.63 $
#
##############################################################################


use Subs::SwiftSub;
use Subs::NominalPointing;
use Util::AttTool;

use Scalar::Util;
use Util::Email;
# Math.pm is in the Util directory (copied from heasoft),
# but isn't part of the Util package.
use Util::Math;

@ISA = ('Subs::NominalPointing', 'Subs::SwiftSub');
use strict;

use Util::SwiftTags;

sub new {
    my $proto=shift;
    my $self=$proto->SUPER::new();
    my $procpar =$self->procpar();

    $self->{DESCRIPTION}="Extracting attitude data";

    return $self;
}

##################
# METHODS:
##################

sub body {
    my $self=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    ##############################################
    # extract attitude files from various sources
    ##############################################
    my @acs = $self->extract_acs("acs_att_tmp");
    my @bat = $filename->get('attitude', 'bat', '', '*');
    my @xrt;
    my @batlc;
    unless( @acs || @bat ){
      $log->entry("No primary attitude data, using secondary sources.");
      push @xrt, $self->extract_xrt("xrt_att_tmp");
      push @bat, $filename->get('attlpd', 'bat', '', '*');
    }else{
      unlink $filename->get('attlpd', 'bat', '', '*');
      my $seq = $jobpar->read('sequence');
      if( $seq%1000==991 ){
        $log->entry("Using xrt attitude data in 991 sequence.");
        push @xrt, $self->extract_xrt("xrt_att_tmp");
      }
    }


    #############################################################
    # merge just the attitude files from the ACS type packets
    ###########################################################
    my $acslist = Util::FITSlist->new(@acs, @bat);
    my @acs_merged;
    unless($acslist->count() == 0) {
      $log->entry("Merging ". join(" ", $acslist->files()) );
      push @acs_merged, 'acs_merged.tmp';
      my $merge_out = $acslist->merge($acs_merged[0]);
      if($merge_out ne $acs_merged[0]) {
        ############################
        # there was only one file
        ###########################
        rename $merge_out, $acs_merged[0];
      }else{
        ##########################
        # delete the files
        ##########################
        foreach my $file ($acslist->files() ) {
          $log->entry("Deleting $file");
          unlink $file;
        }
      }
    } else {
      my $index = 0;
      foreach my $lcfile ($filename->get('batlcatt', 'bat', '', '*'),
              $filename->get('tdrsslcatt', 'bat', '', '*')){

        push @batlc, "bat_lcatt_${index}.tmp";
        Util::FITSfile->new($lcfile, 'ATTITUDE', '[col TIME; QPARAM; POINTING; SOURCE(B)=4]')
                  ->copy($batlc[-1]);
        unlink $lcfile;
        $index++;
      }
    }

    my $attitude = $filename->get('attitude', 's');
    my $tstart = 0;

    # check for attitude override;
    my $att_force = undef;
    {
        my $attdir = $procpar->read('att_force_dir');
        my $seq = $jobpar->read('sequence');
        my $path = "${attdir}/sw${seq}sat.fits";
        if (-f $path) {
            $att_force = $path;
            $log->entry("Found attitude force file: $path");
        }
    }

    #######################################
    # merge all the attitude data
    #######################################
    my $attlist = Util::FITSlist->new(@acs_merged, @batlc, @xrt);
    if($attlist->count() ==0) {
        $log->error(2, 'No attitude files produced, setting ra,dec to 0,0');
        $jobpar->set({ra => 0.0, dec => 0.0});
        return if not $att_force;
    } else {
        $log->entry("Merging ". join(" ", $attlist->files()) );

        my $merged = $attlist->merge($attitude, '', 'TIME', 'QPARAM', 'POINTING', 'SOURCE');

        ###############################################################################
        # Make a seperate extension with just the extra columns from the ACS packets
        ###############################################################################
        if($merged ne $attitude) {
            ############################
            # there was only one file
            ###########################
            rename $merged, $attitude;
            my $inter_file = Util::FITSfile->new($attitude, 'ATTITUDE');
            if($inter_file->find_column('BUS_V')){
              $inter_file->specs('[col TIME; POSITION; FLAGS; BUS_V; SOURCE]');
              $inter_file->append_to($attitude);

              my $delcol = Util::Ftool->new('fdelcol')
                                ->params({infile => $attitude.'[1]',
                                          confirm => 'no',
                                          proceed => 'yes'});

              $delcol->params({colname => 'POSITION'})->run();
              $delcol->params({colname => 'FLAGS'})->run();
              $delcol->params({colname => 'BUS_V'})->run();
            }
        }else{
          if(@acs_merged){
            my $inter_file = Util::FITSfile->new($acs_merged[0], 'ATTITUDE');
            if($inter_file->find_column('BUS_V')){
              $inter_file->specs('[col TIME; POSITION; FLAGS; BUS_V; SOURCE]');
              $inter_file->append_to($attitude);
            }
          }
          ##########################
          # delete the files
          ##########################
          foreach my $file ($attlist->files() ) {
            $log->entry("Deleting $file");
            unlink $file;
          }
        }

        $tstart = $self->normalizeAttitude($attitude);
    }

    if ($att_force) {
        # override attitude- append raw attitude to this file
        my $tmpfile = 'attitude.force';
        Util::FITSfile->new($att_force)
                  ->copy($tmpfile);
        $tstart = $self->normalizeAttitude($tmpfile);

        if (-f $attitude) {
            $log->entry('Appending raw attitude to forced attitude file');
            Util::FITSfile->new($attitude, 'ATTITUDE', q([col EXTNAME='RAW_ATTITUDE']))
                  ->append_to($tmpfile);
            Util::FITSfile->new($attitude, 'ACS_DATA', q([col EXTNAME='RAW_ACS_DATA']))
                  ->append_to($tmpfile);
        }

        rename($tmpfile, $attitude);

        # An E2 error is introduced to ensure operator
        # intervention for the FFA processing
        $log->error(2, "Overrode attitude with $att_force");
    }

    return if not $tstart;

    my $fits = Util::FITSfile->new($attitude, 1);

    my $align_fits = Util::FITSfile->new($filename->fetch_cal('alignment'), 0);
    $fits->keyword('ALGN_NAM', $align_fits->keyword('FILENAME'),
            'Name of alignment teldef file used.');

    #################################################
    # Make settling and pointing gtis
    #################################################
    my $settling = $filename->get('gti', 's', 'st', 0);
    my $pointing = $filename->get('gti', 's', 'po', 0);
    my $not_pointing = $filename->get('gti', 's', 'np', 0);
    my $ontarget = $filename->get('gti', 's', 'ot', 0);
    my $nfis = $filename->get('gti', 's', 'nf', 0);
    my $target_pointing = $filename->get('gti', 's', 'tp', 0);

    if( (grep /^ACS_DATA$/, $fits->list_hdus()) ){
      my $tempfile = 'acs.tmp';
      Util::HEAdas->new('ftsort')
              ->params({infile => $attitude.'[ACS_DATA]',
                outfile => $tempfile,
                unique => 'yes',
                columns => 'TIME, FLAGS'})
          ->run();

      my $maketime = Util::Ftool->new('maketime')
                            ->params({infile => $tempfile.'[ACS_DATA]',
                             compact => 'no'});

      $log->entry("Producing GTIs for settling and pointing.");

      $maketime->params({outfile => $settling,
                         prefr => 0,
                         postfr => 1,
                         expr => 'FLAGS == b10x0xxxx'})
               ->run();

      $maketime->params({outfile => $pointing,
                         prefr => 0,
                         postfr => 0,
                         expr => 'FLAGS == b11x0xxxx'})
               ->run();

      $maketime->params({outfile => $not_pointing,
                         prefr => 1,
                         postfr => 1,
                         expr => 'FLAGS != b11xxxxxx || FLAGS == bxxx1xxxx'})
               ->run();

      if( -s $nfis ){
        Util::Ftool->new('mgtime')
               ->params({ingtis => $pointing .' '. $nfis,
                         outgti => $ontarget,
                         merge => 'AND'})
               ->run();
      }else{
        Util::FITSfile->new($pointing)
                ->copy($ontarget);
      }

      unlink $tempfile if -f $tempfile;
    }

    foreach my $gtifile ($settling, $pointing, $ontarget, $nfis){
      next unless -s $gtifile;

      my $gtifits = Util::FITSfile->new($gtifile, 1);
      unless( $gtifits->nrows() ){
        $log->error([ 1, ATT_NO_GTI ], "GTI file $gtifile is empty. Deleting.");
        unlink $gtifile;
        next;
      }

      $gtifits->keyword('EXTNAME', 'GTI');
      $gtifits->keyword('TSTART', ($gtifits->cols('START')->table())[0] );
      $gtifits->keyword('TSTOP', ($gtifits->cols('STOP')->table())[-1] );
    }

    ##########################################
    # determine the nominal pointing
    # using method inherited from superclass
    ##########################################
    unless( -f $ontarget ){
      if( -f $pointing ) {
        $ontarget = $pointing;
      }else{
        $ontarget = $self->no_gap_gtis($attitude, 100);
      }
    }

    # check if the nominal pointing database has an entry for the sequence
    my $pnt_db = $self->checkDatabaseForNominalPointing;

    if ($pnt_db) {

        $self->storeNominalPointing($pnt_db);

    }
    else {
        my $pnt_sacs = $self->determineRequestedPointing;

        my $pnt_aspect = $self->determineAspect($ontarget);

        if ($pnt_sacs and $pnt_aspect) {

            my $sacs_str = formatPNT($pnt_sacs);
            my $aspect_str = formatPNT($pnt_aspect);

            my $details = "SACSREQ $sacs_str and aspect $aspect_str";
            $log->entry("comparing $details");

            my $ten_arcmin_deg = sprintf('%.2f', 10 / 60);

            if (separation_deg($pnt_sacs, $pnt_aspect) > $ten_arcmin_deg) {
                my $text = "$details positions are inconsistent by more than $ten_arcmin_deg deg";
                my $sequence = $jobpar->read('sequence');
                my $seqprocnum = $jobpar->read('seqprocnum');
                $log->error(2, $text);

                # (NB: Can't include obsdate because that parameter isn't
                # set until a later module)
                $text = $text . "\n\nobject= " . $jobpar->read('object') .
                    "\nsequence= $sequence" .
                    "\ntarget= " . $jobpar->read('target') .
                    "\nprocessing number= ${seqprocnum}" .
                    "\njob title= " . $jobpar->read('job_title') . "\n" ;
                Util::Email::sendEmail(SUB => $self,
                       SUBJECT => "aspect and SACSREQ inconsistent in ${sequence}.${seqprocnum}",
                       TEXT => $text);
            }

            my $rolls = $pnt_sacs->{ROLLS_deg} || [ ];

            my $roll_a = $pnt_aspect->{ROLL_deg};
            my $matched = undef;
            my $nonzero = 0;

            my $max_roll_deg = 2;

            foreach my $roll (@$rolls) {

                if ($roll) {
                    ++$nonzero;
                }

                my $delta_d = deltaRoll_d($roll, $roll_a);
                if ($delta_d < 0) {
                    $log->error(1, "bad input: deltaRoll_d($roll, $roll_a)");
                }
                elsif ($delta_d < $max_roll_deg) {
                    $matched = 1;
                }
            }

            if (not $nonzero) {
                $log->error(1, "no non-zero SACSREQROLL");
            }

            if (not $matched) {
                $log->error(1, "aspect roll not within $max_roll_deg deg of SACSREQROLL");
            }
        }
        elsif ($pnt_aspect) {
            if ($pnt_aspect->{ERROR}) {
                $log->error(2, 'Unable to determine nominal pointing');

                # (NB: Can't include obsdate because that parameter isn't
                # set until a later module)
                my $sequence   = $jobpar->read('sequence');
                my $seqprocnum = $jobpar->read('seqprocnum');
                my $text = "Unable to determine nominal pointing by either aspect or SAQSREQ" .
                    "\n\nobject= " . $jobpar->read('object') .
                    "\nsequence= $sequence" .
                    "\ntarget= " . $jobpar->read('target') .
                    "\nprocessing number= ${seqprocnum}" .
                    "\njob title= " . $jobpar->read('job_title') . "\n" ;
                Util::Email::sendEmail(SUB => $self,
                       SUBJECT => "Unable to determine nominal pointing in ${sequence}.${seqprocnum}",
                       TEXT => $text);

            }
        }
    }

    if($ontarget =~ /\.tmp$/ ) {
        ####################################
        # delete no-gap temporary GTI file
        ####################################
        unlink $ontarget;
    }


    #################################################################
    # GTI of pointings where we are actually pointing in the nominal
    # direction.  Call this 'tp' for 'target pointing.
    #################################################################
    if( -f $pointing ){
      my $maketime = Util::Ftool->new('maketime')
                            ->params({infile => $attitude.'[ATTITUDE]',
                                      compact => 'no',
                                      outfile => 'gti.tmp',
                                      prefr => 1,
                                      postfr => 1,
                                      expr => 'near(POINTING[1],' .$jobpar->read('ra') .',1.0) ' .
                                          ' && near(POINTING[2],' .$jobpar->read('dec') .',1.0)' })
                            ->run();


      Util::Ftool->new('mgtime')
               ->params({ingtis => $pointing .' gti.tmp',
                         outgti => $target_pointing,
                         merge => 'AND'})
               ->run();

      unlink 'gti.tmp';

      my $tpfits = Util::FITSfile->new($target_pointing, 'STDGTI');
      $tpfits->keyword('EXTNAME', 'GTI');
    }

    if (not $att_force) {
      #################################################################
      # Only do aberration correction before velocity adding was turned on
      #################################################################
      my $velocity_add = Util::Date->new('2005-01-31T21:57:00');
      $self->correct_aberration()
        if $tstart < $velocity_add->seconds();
    }

    my $att_proc = $filename->get('attcorr', 'p');
    my $att_corr = Util::HEAdas->new('attjumpcorr')->is_script(1);

    $att_corr->params({infile => $attitude,
               outfile => $att_proc,
               chatter => 3});
    $att_corr->run();

    my $attstatus;
    if( $att_corr->had_error() ){
      $attstatus = 120;
      unlink $att_proc;
    }else{
      $attstatus = 110;
      my $diff = Util::HEAdas->new('ftdiff')
                         ->params({infile1 => $attitude,
                       infile2 => $att_proc});
      $diff->seriousness(0);
      $diff->run();

      my $ndiffs = $diff->parfile()->read('numdiffs');
      if($ndiffs==0){
        $attstatus = 100;
        unlink $att_proc;
      }
    }

    $fits->ext(0);
    $fits->keyword('ATTSTATU', $attstatus, 'Status of corrected attitude files');
    for(my $i=0; $i<3; $i++){
      $fits->ext($i);
      $fits->keyword('ATTFLAG', "'100'", 'Attitude origin: 100=sat/spacecraft');
    }

    if( -f $att_proc){
      my $pattfits = Util::FITSfile->new($att_proc);
      for(my $i=0; $i<3; $i++){
        $pattfits->ext($i);
        $pattfits->keyword('ATTFLAG', "'110'", 'Attitude origin: 110=pat/jump corrected');
      }
    }

    ############################################################
    # calculate attitude/orbit calculated filtering quantities
    # 2016-11-02 moved after attjumpcorr so prefilter can use pat (if available)
    ############################################################
    $self->calculate_att_orb();

} # end of body method



sub normalizeAttitude
{
    my ($self, $attitude) = @_;

    my $log = $self->log;

    ##########################################
    # TIME sort the merged attitude file and
    # remove overdetermined values
    ##########################################
    $log->entry("Normalizing $attitude");

    my $now = Util::Date->new();
    my $Tnow = $now->date() .'T'. $now->time();
    my $fits = Util::FITSfile->new($attitude, 0);
    $fits->keyword('DATE', $Tnow);
    $fits->keyword('MJDREFI', 51910, 'MJD reference day');
    $fits->keyword('MJDREFF', 7.428703700000000E-04, 'MJD reference (fraction of day)');

    $fits->ext(1);
    $fits->sort();

    if($fits->nhdus() > 2){
      $fits->ext(2);
      $fits->cols('TIME');
      $fits->sort('unique')
           ->keyword('EXTNAME', 'ACS_DATA');
    }

    Util::AttTool->new("att_thin")
             ->command_line($attitude)
             ->run();

    $log->entry("normalizing quaternions");
    my $colfilter = '[col *; QPARAM = QPARAM / sqrt(sum(QPARAM*QPARAM))]';
    Util::HEAdas->new('ftcopy')
            ->params({
              infile => $attitude . $colfilter,
              outfile => 'attitude.tmp',
             })
        ->run;

    ###############################################################
    # Run diff so that we have a record in the log of what changed
    ###############################################################
    my $diff = Util::HEAdas->new('ftdiff')
                       ->params({infile1 => "$attitude\[1]",
                     infile2 => "attitude.tmp\[1]",
                     reltol => 10e-8});

    $diff->run();

    rename('attitude.tmp', $attitude);

    ##########################################
    # set TSTART, TSTOP, and UTCFINIT in the merged file
    ##########################################
    $fits = Util::FITSfile->new($attitude, 1)
                          ->cols("TIME");

    my $nrows = $fits->nrows();
    my $tstart = $fits->rows(1     )->table();
    my $start_date = Util::Date->new($tstart);
    my $start_str = $start_date->date() . 'T' . $start_date->time();
    my $tstop  = $fits->rows($nrows)->table();
    my $stop_date = Util::Date->new($tstop);
    my $stop_str = $stop_date->date() . 'T' . $stop_date->time();
    my $utcfinit = undef;

    {
        # determine UTCFINIT at TSTART
        my $filename = $self->filename;
        my $time_file = $filename->get('timedata', 'swift', '', 0);
        my (@utcf_times, @utcf);
        if (-f $time_file) {
            my $time_fits = Util::FITSfile->new($time_file, 'UTCF');
            @utcf_times = $time_fits->cols('TIME')->table();
            @utcf = $time_fits->cols('UTCF')->table();
        }

        if (@utcf) {
            my $itime = 0;
            while ($tstart > $utcf_times[$itime] &&
                    $itime < $#utcf_times) { $itime++ }
            $utcfinit = $utcf[$itime];
        }

        if (defined($utcfinit)) {
            $log->entry("found UTCFINIT=$utcfinit at $tstart");
        }
        else {
            $log->error(1, "unknown UTCFINIT at $tstart");
        }
    }

    my $nhdus = $fits->nhdus;
    for (my $i = 0; $i < $nhdus; ++$i) {

        $fits->ext($i);

        $fits->keyword('DATE', $Tnow);
        $fits->keyword('MJDREFI', 51910, 'MJD reference day');
        $fits->keyword('MJDREFF', 7.428703700000000E-04, 'MJD reference (fraction of day)');

        $fits->keyword('TSTART', $tstart, 'Time of first attitude record');
        $fits->keyword('TSTOP' , $tstop, 'Time of last attitude record');
        $fits->keyword('DATE-OBS', $start_str, 'Date corresponding to TSTART');
        $fits->keyword('DATE-END', $stop_str, 'Date corresponding to TSTOP');

        if (defined($utcfinit)) {
            $fits->keyword('UTCFINIT', $utcfinit, '[s] UTCF at TSTART');
        }
    }

    return $tstart;
}


#############################################################################
# Extract an attitude file from the ACS packets
#############################################################################
sub extract_acs {
    my $self=shift;
    my $base=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    my @list=();


    ##############################
    # get the ACS packets file
    ##############################
    my @acs = ( $filename->get("telemetry", "*", "acs", "*"),
            $filename->get("telemetry", "*", "head2", 485) );
    unless(@acs && -f $acs[0]) {
        $log->entry("No ACS packets available");
        return @list;
    }


    ##########################
    # run acs2fits
    ##########################
    my $index=0;
    foreach my $ccsds_file (@acs) {
        #######################
        # get the file name
        #######################
        my $file = "$base.$index";
        $log->entry("Creating $file from $ccsds_file");
        unlink $file;

        #############################
        # set up the extraction tool
        #############################
        my $teldef = $filename->fetch_cal('alignment');
        my $acs2fits=Util::AttTool->new('acs2fits');

        $acs2fits->command_line("-infile $ccsds_file",
                "-alignfile $teldef",
                                "-outfile $file");

        #####################################
        # run the tool and check for errors
        #####################################
        $acs2fits->run();
        if($acs2fits->had_error() ) {
            #############
            # error
            #############
            if(-f $file ) {
                $log->entry("Removing $file");
                unlink $file;
            }
        } else {
            #############
            # success
            #############
            push @list, ($file);
        }


    } # end of loop over input files

    return @list;

} # end of extract_acs method


#############################################################################
# Extract an attitude file from the XRT FITS files
#############################################################################
sub extract_xrt {
    my $self=shift;
    my $base=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    my @list=();

    ###########################################################
    # types of XRT HK files with attitude data in them
    ###########################################################
    my @modes = ('48a', '4e0', '4e1', '4e2', '4f0', '500', '534', '536', '53a', '540');
###    "48Ah", "4E0h", "4E1h", "4E2h",
###    "4F0h", "500h", "534h", "536h", "53Ah", "540h");

    my @files=();
    my @head_files = ( $filename->get('hk', 'xrt', 'hd', '*') );
    push @files, map $_ .= "[1]", @head_files;

    my $eng_hk = $filename->get('enhk', 'xrt', '', 0);
    if( -f $eng_hk ){
      my $fitsEng = Util::FITSfile->new($eng_hk);
      my $nhdus = $fitsEng->nhdus();
      for(my $iext=1; $iext < $nhdus; $iext++){
    $fitsEng->ext($iext);
    my $name = $fitsEng->keyword('EXTNAME');
    if( grep(/^hk${name}x$/, @modes) ){
      push @files, $eng_hk . "[$name]";
        }
      }
    }


    unless(@files) {
        $log->entry("No attitude data in the XRT FITS files");
        return @list;
    }


    ###################################
    # extract the attitude files
    ###################################
    my $index=0;
    foreach my $xrt_file (@files) {
        #######################
        # get the file name
        #######################
        my $file = "$base.".$index++;
        $log->entry("Creating $file from $xrt_file");
        unlink $file;

        #############################
        # set up the extraction tool
        #############################
    my $teldef = $filename->fetch_cal('alignment');
        my $tool=Util::AttTool->new('xrtatt');

        $tool->command_line("-infile '$xrt_file'",
                "-alignfile $teldef",
                            "-outfile $file");

        #####################################
        # run the tool and check for errors
        #####################################
        $tool->run();
        if($tool->had_error() ) {
            #############
            # error
            #############
            if(-f $file ) {
                $log->entry("Removing $file");
                unlink $file;
            }
        } else {
            #############
            # success
            #############
            push @list, ($file);
        }


    } # end of loop over input files

    return @list;


} # end of extract_xrt method



#############################################################################
# create a GTI file from an attitude file which excludes gaps in the
# attitude records
#############################################################################
sub no_gap_gtis {
    my $self = shift;
    my $att = shift;
    my $gap = shift;

    my $log     =$self->log();
    my $filename=$self->filename();

    $log->entry("Creating GTI file which excludes attitude ".
                "record gaps > $gap s");

    #############################################
    # read the time records in the attitude file
    #############################################
    my $fits = Util::FITSfile->new($att, 1);
    my @time = $fits->cols("TIME")->table();

    unless(@time) {
        $log->entry("No attitude records");
        return "";
    }

    #######################################
    # open an ASCII file to hold the data
    #######################################
    my $data = "att_gtis.tmp";
    unlink $data;
    open DATA, ">$data";


    my $start=$time[0];
    my $i;
    my $nrows = $fits->nrows();
    for($i=0; $i<$nrows; $i++) {

        $start = $time[$i];

        ###############################################################
        # skip over consecutive rows to find the start of the next gap
        ###############################################################
        while ($i < $nrows && $time[$i] - $time[$i-1] < $gap ) { $i++ }

        ####################
        # record the GTI
        ####################
        print DATA "$start $time[$i-1]\n";

        ###############################################################
        # skip over the gap
        ###############################################################
        while ($i < $nrows && $time[$i] - $time[$i-1] >= $gap ) { $i++ }


    } # end of loop over attitude records

    close DATA;


    my $header = "att_gtis_header.tmp";
    open HEADER, ">$header";
    print HEADER "START 1D\n";
    print HEADER "STOP  1D\n";
    close HEADER;

    # my $gti = "att_gtis_fits.tmp";
    my $gti = $filename->get('gti', 's', 'ng', 0);

    my $fcreate = Util::Ftool->new("fcreate");

    $fcreate->params({cdfile   => $header,
                      datafile => $data,
                      outfile  => $gti,
                      headfile => " ",
                      tbltype  => "binary",
                      nskip    => 0,
                      nrows    => 0,
                      history  => "yes",
                      morehdr  => 0,
                      extname  => "GTI",
                      anull    => " ",
                      inull    => 0,
                      clobber  => "yes"})
            ->run();

    unlink $data;
    unlink $header;

    unless($fcreate->had_error() ) { return $gti }
    else                           { return "" }

} # end of no_gap_gtis;

#############################################################################
# Calculate attitude/orbit derived filter quantities
#############################################################################
sub calculate_att_orb {
    my $self=shift;
    my @columns=@_;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    $log->entry("Calculating attitude/orbit filtering quantities");

    #######################################
    # make sure we have an attitude file
    # othewise there's no sense bothering
    #######################################
    my $attitude = $filename->get('attcorr', 'p');
    unless (-f $attitude) {
        $attitude = $filename->get('attitude', 's');
        unless (-f $attitude) {
            $log->entry("No attitude file available");
            return;
        }
    }

    $log->entry("Running prefilter with $attitude attitude");

    ##################################################################
    # assemble the mission zero time into a string to feed prefilter
    ##################################################################
    my $epoch = $procpar->read("refdate")."T".$procpar->read("reftime");

    my $fits = Util::FITSfile->new($attitude, 1);
    ##############################################
    # Add a buffer of 60 seconds to start and stop
    ##############################################
    my $tstart = $fits->keyword("TSTART") - 60.0;
    my $tstop  = $fits->keyword("TSTOP") + 60.0;

    my $attorb = $filename->get('attorb', 's');
    unlink $attorb;

    my $cols = 'ALL';
    $cols = join(' ', @columns) if @columns;

    my $orbname = $filename->fetch_orbit();

    #####################################
    # set up and run prefilter
    #####################################
    my $prefilter = Util::HEAdas->new("prefilter");
    $prefilter->params({outname   => $attorb,
                        columns   => $cols,
                        orbmode   => "TLE_TEXT2",
                        orbname   => $orbname,
                        attname   => $attitude,
                        alignfile => 'CALDB',
                        leapname  => $filename->fetch_cal("leapsec"),
                        rigname   => $filename->fetch_cal("rigidity"),
                        start     => $tstart,
                        end       => $tstop,
                        interval  => 1.0, # was 30.0
                        attextrap => '32',
                        origin    => 'GSFC',
                        ranom     => $jobpar->read("ra"),
                        decnom    => $jobpar->read("dec"),
                        missepoch => $epoch,
                        compcols  => "TIME",
                        compapplyquat => "no",
                        chatter   => 4,
                        clobber   => "yes" });

    $log->entry("Running ". $prefilter->name() );
    $prefilter->run();



} # end of calculate_att_orb

#############################################################################
# Calculate attitude/orbit derived filter quantities
#############################################################################
sub correct_aberration {
    my $self=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    $log->entry('Doing aberration correction of attitude file');

    #######################################
    # make sure we have an attitude file
    # othewise there's no sense bothering
    #######################################
    my $attitude = $filename->get('attitude', 's');
    unless( -f $attitude ) {
        $log->entry('No attitude data available');
        return;
    }

    $log->entry('Preliminary prefilter run');
    $self->calculate_att_orb('TIME', 'POSITION', 'VELOCITY');

    my $aberration = Util::HEAdas->new('aberrator')
                                 ->params({infile => $attitude,
                       orbfile => $filename->get('attorb', 's'),
                       alignfile => 'CALDB'
                      });

    $log->entry('Run aberrator');
    $aberration->run();

} # end of correct_aberration



sub determineAspect
{
    my ($self, $ontarget) = @_;

    my $problem = 0;
    my $level = 1;

        # determine_pointing is found in comN.N/Subs/NominalPointing.pm.
        # It stores the aspect in the ra, dec, and roll parameters of the
        # jobpar file.
    $self->determine_pointing($ontarget, ERROR_LEVEL => $level,
                  PROBLEM_REF => \$problem);

    my $jobpar = $self->jobpar;
    my $ra = $jobpar->read('ra');
    my $dec = $jobpar->read('dec');
    my $roll = $jobpar->read('roll');

    my %aspect = (
        TAG => 'aspect',
        RA_deg => $ra,
        DEC_deg => $dec,
        ROLL_deg => $roll,
        ERROR => $problem,
    );

    return \%aspect;
}


sub checkDatabaseForNominalPointing
{
    my ($self) = @_;

    my $procpar = $self->procpar;
    my $jobpar = $self->jobpar;
    my $log = $self->log;

    my $sequence = $jobpar->read('sequence');

    my $pnt = undef;

    my $rdbbin = $ENV{RDB_BIN} || '/aps/tools/rdb';
    my $rdbtables = $ENV{RDB_TABLES_DIR} || '/aps/db/rdb/tables';

    my $rdbpath = $rdbtables . '/sw/nominal_pnt.rdb';
    if (not -f $rdbpath) {
        $log->error(1, "missing nominal pointing table");
    }
    else {
        my @out = qx($rdbbin/row < $rdbpath sequence eq $sequence | $rdbbin/headchg -del);
        if (@out == 1) {

        my $record = $out[0];
        chomp($record);

        $log->entry("nominal pointing override: $record");

        my ($sequence, $ra, $dec, $roll) = split(/\t/, $record);

        my $ok = 0;
        foreach my $x ($sequence, $ra, $dec, $roll) {
            if (Scalar::Util::looks_like_number($x)) {
            ++$ok;
            }
            else {
            $log->error(1, "bad value '$x' for $sequence");
            }
        }

        if ($ok == 4) {
            $pnt = {
            TAG => 'RDB',
            RA_deg => $ra,
            DEC_deg => $dec,
            ROLL_deg => $roll,
            };
        }
        }
        else {
        $log->entry("no nominal pointing in DB for $sequence");
        }
    }

    return $pnt;
}


sub formatPNT
{
    my ($p) = @_;

    my $roll;
    if (ref($p->{ROLLS_deg})) {
        my @rolls = map { sprintf('%.1f', $_) } @{ $p->{ROLLS_deg} };
        $roll = join('/', @rolls);
    }
    else {
        $roll = sprintf('%.1f', $p->{ROLL_deg});
    }
    my $str = sprintf('RA=%.3f Dec=%.3f Roll=%s [deg]',
            $p->{RA_deg}, $p->{DEC_deg}, $roll);

    return $str;
}


sub d2r
{
    my ($deg) = @_;
    my $rad = Math::toRadians($deg);
    return $rad;
}


sub separation_deg
{
    my ($p1, $p2) = @_;
    my $u1 = Math::rd2unit(d2r($p1->{RA_deg}), d2r($p1->{DEC_deg}));
    my $u2 = Math::rd2unit(d2r($p2->{RA_deg}), d2r($p2->{DEC_deg}));

    my $rad = Math::u3angle($u1, $u2);
    my $deg = Math::toDegrees($rad);
    return $deg;
}


sub storeNominalPointing
{
    my ($self, $pnt) = @_;

    my $log = $self->log;
    my $jobpar = $self->jobpar;

    my $ra = $pnt->{RA_deg} || 0;
    my $dec = $pnt->{DEC_deg} || 0;
    my $roll = $pnt->{ROLL_deg} || 0;

    my ($euler1, $euler2, $euler3)
         = convertRADecRollToSwiftEuler($ra, $dec, $roll);

    #############################################
    # convert to galactic coordinates
    #############################################
    my $ftcoco = Util::HEAdas->new('ftcoco')->is_script(1)->params({
        infile   => 'NONE',
        outfile  => '',
        incoord  => 'R',
        outcoord => 'G',
        lon1     => $ra,
        lat1     => $dec,
        radecsys => 'FK5',
        chatter  => 0,
        clobber  => 'no',
        history  => 'no'
        })->run();

    # Get the output Lii/Bii params if no error
    my $glon = 0.0;
    my $glat = 0.0;
    if ($ftcoco->had_error()) {
        $log->error( 2, "Error in ftcoco!" );
    }
    else {
        $glon = $ftcoco->parfile()->read('lon2', 1e-6);
        $glat = $ftcoco->parfile()->read('lat2', 1e-6);
    }

    ####################################################
    # log the results and write them to the jobpar file
    ####################################################
    $log->entry("Nominal spacecraft Euler angles: ".
        "Phi=$euler1 Theta=$euler2 Psi=$euler3");

    $log->entry("Nominal telescope pointing: ".
        "R.A.=$ra Dec=$dec Roll Angle=$roll");

    $log->entry("Galactic Coordinates: ".
        "Lii=$glon Bii=$glat");

    $jobpar->set({
        ra     => $ra,
        dec    => $dec,
        roll   => $roll,
        euler1 => $euler1,
        euler2 => $euler2,
        euler3 => $euler3,
        glon   => $glon,
        glat   => $glat });
}


sub determineRequestedPointing
{
    my ($self) = @_;

    my $log      = $self->log;
    my $filename = $self->filename;
    my $procpar  = $self->procpar;
    my $jobpar   = $self->jobpar;

    my $senhk = $filename->get('scenhk', 'swift', '', 0);
    if (not -f $senhk) {
        $log->error(1, "unable to determine SACSREQ without engineering housekeeping");
        return;
    }

    my $outfile = 'sacsreq.fits';

    my $obsid = $jobpar->read('sequence');
    Util::FITSfile->new($senhk, 0)->keyword(OBS_ID => $obsid);

    my $swacsreq = Util::HEAdas->new('swacsreq')
        ->is_script(1)
        ->params({
        infile   => $senhk,
        outfile  => $outfile,
        chatter  => 3,
        clobber  => 'yes',
        history  => 'no'
        })->run();

    if (not -f $outfile) {
        $log->error(1, "unable to determine requested pointing from $senhk");
        return;
    }

    my %keys = Util::FITSfile->new($outfile, 0)->keywords;

    my ($ra, $dec);
    my @roll;
    foreach my $key (keys(%keys)) {

        next if ($key !~ /^SW/);

        my $value = $keys{$key};

        if ($key eq 'SWRA') {
            $ra = $value;
        }
        elsif ($key eq 'SWDEC') {
            $dec = $value;
        }
        elsif ($key =~ /^SWROLL/) {
            push(@roll, $value);
        }
    }

    my $sacsreq = undef;

    if (defined($ra) and defined($dec) and @roll) {
        $sacsreq = {
            RA_deg => $ra,
            DEC_deg => $dec,
            ROLL_deg => $roll[0],
            ROLLS_deg => \@roll,
        };
        $log->entry("determined requested pointing: RA=$ra Dec=$dec roll=$roll[0]");
    }
    else {
        $log->error(1, "unable to determine requested pointing from $senhk: missing RA/Dec/roll keys");
    }

    return $sacsreq;
}


sub testCelestialToEuler
{
    testOne(325.4736, 2.669, 62.4474);
    # => 56.8655, 62.47986, -93.01

    testOne(299.5535, 35.1857, 307.5312);
    # => -174.3233, 60.1397, 131.639

    testOne(186.44, 16.321, 293.61)
    # => 89.4374, 67.3954, 107.722
}


sub testOne
{
    my @p = @_;

    print join("\n\t", 'pointing:', @p), "\n";

    my @e = convertRADecRollToSwiftEuler(@p);
    print join("\n\t", 'euler:', @e), "\n";
}


sub convertRADecRollToSwiftEuler
{
    my ($ra_d, $dec_d, $roll_d) = @_;

    $roll_d *= -1;

    my $d2r = Math::toRadians(1);
    my $r2d = 1 / $d2r;

    my @e_d = ($ra_d, 90 - $dec_d, $roll_d + 90);
    my @e_r = map { $_ * $d2r } @e_d;

    my $cos_phi = cos($e_r[0]);
    my $sin_phi = sin($e_r[0]);

    my $cos_theta = cos($e_r[1]);
    my $sin_theta = sin($e_r[1]);

    my $cos_psi = cos($e_r[2]);
    my $sin_psi = sin($e_r[2]);

    # celestial unit vectors
    my @celux = (
        $cos_phi*$cos_theta*$cos_psi - $sin_phi*$sin_psi,
        $sin_phi*$cos_theta*$cos_psi + $cos_phi*$sin_psi,
        -$sin_theta*$cos_psi
    );

    my @celuy = (
        -$cos_phi*$cos_theta*$sin_psi - $sin_phi*$cos_psi,
        -$sin_phi*$cos_theta*$sin_psi + $cos_phi*$cos_psi,
        $sin_theta*$sin_psi
    );

    my @celuz = (
        $cos_phi*$sin_theta,
        $sin_phi*$sin_theta,
        $cos_theta,
    );

    # spacecraft unit vectors
    my @scx = @celuz;
    my @scy = @celux;
    my @scz = @celuy;

    # spacecraft euler angles
    my $scz2 = $scz[2];
    # range check and circumvent exception
    if (abs($scz2) > 1) {
        if (abs($scz2) - 1 > 1e-6) {
            print "invalid unit vector: scz2=$scz2\n";
        }
        $scz2 = $scz2 < 0 ? -1 : 1;
    }

    my @sce_r = (
        POSIX::atan2($scz[1], $scz[0]),
        POSIX::acos($scz2),
        POSIX::atan2($scy[2], -$scx[2])
    );

    my @sce_d = map { $r2d * $_ } @sce_r;

    return @sce_d;
}


sub deltaRoll_d
{
    my ($x, $y) = @_;

    if ($x < 0 || $x > 360) {
        # carp("invalid roll $x");
        return -1;
    }

    if ($y < 0 || $y > 360) {
        # carp("invalid roll $x");
        return -1;
    }

    my $deg = undef;

    my $q = abs($x - $y);
    if ($q > 180) {
        if ($x < $y) {
            $deg = $x + 360 - $y;
        }
        else {
            $deg = $y + 360 - $x;
        }
    }
    else {
        $deg = $q;
    }

    return $deg;
}


1;