diff --git a/starcheck/calc_ccd_temps.py b/starcheck/calc_ccd_temps.py index 10b1005d..778d1a97 100755 --- a/starcheck/calc_ccd_temps.py +++ b/starcheck/calc_ccd_temps.py @@ -42,6 +42,8 @@ import proseco.characteristics as proseco_char from xija.get_model_spec import get_xija_model_spec +from starcheck.utils import replace_with_obc_quats + from starcheck import __version__ as version MSID = {'aca': 'AACCCDPT'} @@ -344,6 +346,8 @@ def get_week_states(rltt, sched_stop, bs_cmds, tlm): # Add in the backstop commands cmds = cmds.add_cmds(bs_cmds) + cmds = replace_with_obc_quats(cmds) + # Get the states for available commands. This automatically gets continuity. state_keys = ['obsid', 'pitch', 'q1', 'q2', 'q3', 'q4', 'eclipse'] states = kadi_states.get_states(cmds=cmds, start=init_tlm_time, stop=sched_stop, diff --git a/starcheck/check_ir_zone.py b/starcheck/check_ir_zone.py index 24083b96..8f92f97b 100644 --- a/starcheck/check_ir_zone.py +++ b/starcheck/check_ir_zone.py @@ -3,6 +3,7 @@ import kadi.commands import kadi.commands.states from cxotime import CxoTime +from starcheck.utils import replace_with_obc_quats def ir_zone_ok(backstop_file, out=None): @@ -19,6 +20,8 @@ def ir_zone_ok(backstop_file, out=None): # Get the states for available commands. This automatically gets continuity. state_keys = ['pcad_mode', 'obsid'] + bs_cmds = replace_with_obc_quats(bs_cmds) + states = kadi.commands.states.get_states(cmds=bs_cmds, start=rltt, stop=sched_stop, state_keys=state_keys, merge_identical=True) diff --git a/starcheck/src/lib/Ska/Parse_CM_File.pm b/starcheck/src/lib/Ska/Parse_CM_File.pm index 5afe6808..ef7f95d8 100644 --- a/starcheck/src/lib/Ska/Parse_CM_File.pm +++ b/starcheck/src/lib/Ska/Parse_CM_File.pm @@ -31,20 +31,6 @@ our @EXPORT = qw(); our @EXPORT_OK = qw( ); %EXPORT_TAGS = (all => \@EXPORT_OK); -############################################################### -sub rel_date2time { -############################################################### - - # Return seconds when suppled a "relative datetime" of the - # format 000:00:00:00.000 (DOY:HH:MM:SS.sss). - my $date = shift; - - # The old code here uses reverse to just ignore a year if - # included in the string. - my ($sec, $min, $hr, $doy) = reverse split ":", $date; - return ($doy * 86400 + $hr * 3600 + $min * 60 + $sec); -} - ############################################################### sub TLR_load_segments { ############################################################### @@ -531,10 +517,6 @@ sub DOT { %{ $dot{$_} } = parse_params($command{$_}); $dot{$_}{time} = date2time($dot{$_}{TIME}) if ($dot{$_}{TIME}); - # MANSTART is in the dot as a "relative" time like "000:00:00:00.000", so just pass it - # to the rel_date2time routine designed to handle that. - $dot{$_}{time} += rel_date2time($dot{$_}{MANSTART}) - if ($dot{$_}{TIME} && $dot{$_}{MANSTART}); $dot{$_}{cmd_identifier} = "$dot{$_}{anon_param1}_$dot{$_}{anon_param2}" if ($dot{$_}{anon_param1} and $dot{$_}{anon_param2}); $dot{$_}{linenum} = $linenum{$_}; @@ -755,134 +737,6 @@ sub PS { return @ps; } -############################################################### -sub MM { - - # Parse maneuver management (?) file -############################################################### - # This accepts a reference to a hash as the only argument - # the return type may be specified in the hash as 'hash' or 'array' - # default return is hash - # With regard to the return data: - - my $arg_ref = shift; - my $mm_file = $arg_ref->{file}; - my $ret_type = 'hash'; - - if (defined $arg_ref->{ret_type}) { - $ret_type = $arg_ref->{ret_type}; - } - - my $manvr_offset = 10; # seconds expected from AONMMODE to AOMANUVR - - my @mm_array; - - my $mm_text = io($mm_file)->slurp; - - # split the file into maneuvers - my @sections = split(/MANEUVER\sDATA\sSUMMARY\n/, $mm_text); - - # ignore pieces of the file without ATTITUDES - my @good_sect = grep { /INITIAL|FINAL/ } @sections; - - my $int_obsid = 'IN_IA'; - for my $entry (@good_sect) { - - # only keep the relevant bits of each entry (before OUTPUT DATA) - my @para = split(/\n\n/, $entry); - my @attitudes = grep { /ATTITUDE/ } @para; - if (scalar(@attitudes) > 2) { - croak("Maneuver Summary has too many attitudes in section\n"); - } - - # where final or initial attitude may be an intermediate attitude - my @output_data_match = grep { /OUTPUT\sDATA/ } @para; - my $output_data = $output_data_match[0]; - my $initial_attitude = $attitudes[0]; - my $final_attitude = $attitudes[1]; - - my %manvr_hash; - $manvr_hash{initial_obsid} = $1 - if ($initial_attitude =~ /INITIAL ID:\s+(\S+)\S\S/); - $manvr_hash{final_obsid} = $1 if ($final_attitude =~ /FINAL ID:\s+(\S+)\S\S/); - $manvr_hash{start_date} = $1 - if ($initial_attitude =~ /TIME\s*\(GMT\):\s+(\S+)/); - $manvr_hash{stop_date} = $1 if ($final_attitude =~ /TIME\s*\(GMT\):\s+(\S+)/); - $manvr_hash{ra} = $1 if ($final_attitude =~ /RA\s*\(deg\):\s+(\S+)/); - $manvr_hash{dec} = $1 if ($final_attitude =~ /DEC\s*\(deg\):\s+(\S+)/); - $manvr_hash{roll} = $1 if ($final_attitude =~ /ROLL\s*\(deg\):\s+(\S+)/); - $manvr_hash{dur} = $1 if ($output_data =~ /Duration\s*\(sec\):\s+(\S+)/); - $manvr_hash{angle} = $1 - if ($output_data =~ /Maneuver Angle\s*\(deg\):\s+(\S+)/); - my @quat = ($1, $2, $3, $4) - if ($final_attitude =~ /Quaternion:\s+(\S+)\s+(\S+)\s+(\S+)\s+(\S+)/); - $manvr_hash{q1} = $quat[0]; - $manvr_hash{q2} = $quat[1]; - $manvr_hash{q3} = $quat[2]; - $manvr_hash{q4} = $quat[3]; - $manvr_hash{tstart} = date2time($manvr_hash{start_date}); - $manvr_hash{tstop} = date2time($manvr_hash{stop_date}); - - # let's just add those 10 seconds to the summary tstart so it lines up with - # AOMANUVR in backstop - $manvr_hash{tstart} += $manvr_offset; - $manvr_hash{start_date} = time2date($manvr_hash{tstart}); - - # clean up obsids (remove prepended 0s) - if (defined $manvr_hash{initial_obsid}) { - $manvr_hash{initial_obsid} =~ s/^0+//; - } - - # use a dummy or the last initial attitude if there isn't one - else { - $manvr_hash{initial_obsid} = $int_obsid; - } - - if (defined $manvr_hash{final_obsid}) { - $manvr_hash{final_obsid} =~ s/^0+//; - } - else { - $int_obsid = $manvr_hash{initial_obsid} . '_IA'; - $manvr_hash{final_obsid} = $int_obsid; - } - - $manvr_hash{obsid} = $manvr_hash{final_obsid}; - push @mm_array, \%manvr_hash; - - } - - # create a manvr_dest key to record the eventual destination of - # all manvrs. - for my $i (0 .. $#mm_array) { - - # by default the destination is just the final_obsid - $mm_array[$i]->{manvr_dest} = $mm_array[$i]->{final_obsid}; - - # but if the final_obsid has the string that indicates it is - # an intermediate attitude, loop through the rest of the manvrs - # until we hit one that isn't an intermediate attitude - next unless ($mm_array[$i]->{final_obsid} =~ /_IA/); - for my $j ($i .. $#mm_array) { - next if ($mm_array[$j]->{final_obsid} =~ /_IA/); - $mm_array[$i]->{manvr_dest} = $mm_array[$j]->{final_obsid}; - last; - } - } - - if ($ret_type eq 'array') { - return @mm_array; - } - - my %mm_hash; - for my $manvr (0 ... $#mm_array) { - my $obsid = $mm_array[$manvr]->{final_obsid}; - $mm_hash{$obsid} = $mm_array[$manvr]; - } - - return %mm_hash; - -} - ##*************************************************************************** sub mechcheck { ##*************************************************************************** diff --git a/starcheck/src/lib/Ska/Starcheck/Obsid.pm b/starcheck/src/lib/Ska/Starcheck/Obsid.pm index b76429fc..742068c6 100644 --- a/starcheck/src/lib/Ska/Starcheck/Obsid.pm +++ b/starcheck/src/lib/Ska/Starcheck/Obsid.pm @@ -284,15 +284,42 @@ sub set_target { ################################################################################## my $self = shift; - my $manvr = find_command($self, "MP_TARGQUAT", -1); # Find LAST TARGQUAT cmd - ($self->{ra}, $self->{dec}, $self->{roll}) = - $manvr - ? quat2radecroll($manvr->{Q1}, $manvr->{Q2}, $manvr->{Q3}, $manvr->{Q4}) - : (undef, undef, undef); + # Get quat from MP_TARGQUAT (backstop) command. + my $c = find_command($self, "MP_TARGQUAT", -1); # Find LAST TARGQUAT cmd + if (defined $c) { - $self->{ra} = defined $self->{ra} ? sprintf("%.6f", $self->{ra}) : undef; - $self->{dec} = defined $self->{dec} ? sprintf("%.6f", $self->{dec}) : undef; - $self->{roll} = defined $self->{roll} ? sprintf("%.6f", $self->{roll}) : undef; + # Check if it has numeric issues for Quaternion package + my $one_minus_xn2 = 1 - (2 * ($c->{Q1} * $c->{Q3} - $c->{Q2} * $c->{Q4}))**2; + if ($one_minus_xn2 < -1e-12){ + push @{ $self->{warn} }, + sprintf("Backstop 4-element quat not normalized for one_minus_xn2 test\n"); + } + # Compute 4th component (as only first 3 are uplinked) and renormalize. + # Intent is to match OBC Target Reference subfunction + my $q4_obc = sqrt(abs(1.0 - $c->{Q1}**2 - $c->{Q2}**2 - $c->{Q3}**2)); + my $norm = sqrt($c->{Q1}**2 + $c->{Q2}**2 + $c->{Q3}**2 + $q4_obc**2); + if (abs(1.0 - $norm) > 1e-6) { + push @{ $self->{warn} }, + sprintf("Uplink quaternion norm value $norm is too far from 1.0\n"); + } + my @c_quat_norm = + ($c->{Q1} / $norm, $c->{Q2} / $norm, $c->{Q3} / $norm, $q4_obc / $norm); + ($self->{ra}, $self->{dec}, $self->{roll}) = quat2radecroll(@c_quat_norm); + ($self->{OBC_Q1}, $self->{OBC_Q2}, $self->{OBC_Q3}, $self->{OBC_Q4}) = + @c_quat_norm; + } + else { + ($self->{ra}, $self->{dec}, $self->{roll}) = (undef, undef, undef); + ($self->{OBC_Q1}, $self->{OBC_Q2}, $self->{OBC_Q3}, $self->{OBC_Q4}) = + (undef, undef, undef, undef); + } + + + # Check abs declination MP guideline + if ((defined $self->{dec}) and (abs($self->{dec}) >= 89.7)){ + push @{ $self->{warn} }, + sprintf("Target abs(Dec) is %.1f degrees > 89.7 limit \n", abs($self->{dec})); + } } @@ -325,6 +352,17 @@ sub find_command { return undef; } +sub first_maneuver_after { + my $man_time = shift; + my $mm = shift; + MANVR: + foreach my $m (@{$mm}) { + next MANVR if $m->{tstart} < $man_time; + return $m; + } + return undef; +} + ################################################################################## sub set_maneuver { # @@ -339,72 +377,58 @@ sub set_maneuver { while ($c = find_command($self, "MP_TARGQUAT", $n++)) { $found = 0; - foreach my $m (@{$mm}) { - my $manvr_obsid = $m->{final_obsid}; - -# where manvr_dest is either the final_obsid of a maneuver or the eventual destination obsid - # of a segmented maneuver - if ( ($manvr_obsid eq $self->{dot_obsid}) - && abs($m->{q1} - $c->{Q1}) < 1e-7 - && abs($m->{q2} - $c->{Q2}) < 1e-7 - && abs($m->{q3} - $c->{Q3}) < 1e-7) + + my $m = first_maneuver_after($c->{time}, $mm); + if ( (defined $m) + && abs($m->{q1} - $c->{Q1}) < 1e-7 + && abs($m->{q2} - $c->{Q2}) < 1e-7 + && abs($m->{q3} - $c->{Q3}) < 1e-7) + { + $found = 1; + foreach my $key ( + qw(initial_obsid final_obsid dur angle ra dec roll q1 q2 q3 q4 tstart start_date tstop stop_date) + ) { - $found = 1; - foreach (keys %{$m}) { - $c->{$_} = $m->{$_}; - } + $c->{$key} = $m->{$key}; + } # Set the default maneuver error (based on WS Davis data) and cap at 85 arcsec - $c->{man_err} = (exists $c->{angle}) ? 35 + $c->{angle} / 2. : 85; - $c->{man_err} = 85 if ($c->{man_err} > 85); - - # Get quat from MP_TARGQUAT (backstop) command. - # Compute 4th component (as only first 3 are uplinked) and renormalize. - # Intent is to match OBC Target Reference subfunction - my $q4_obc = sqrt(abs(1.0 - $c->{Q1}**2 - $c->{Q2}**2 - $c->{Q3}**2)); - my $norm = sqrt($c->{Q1}**2 + $c->{Q2}**2 + $c->{Q3}**2 + $q4_obc**2); - if (abs(1.0 - $norm) > 1e-6) { - push @{ $self->{warn} }, - sprintf( - "Uplink quaternion norm value $norm is too far from 1.0\n"); - } - - my @c_quat_norm = ( - $c->{Q1} / $norm, - $c->{Q2} / $norm, - $c->{Q3} / $norm, - $q4_obc / $norm - ); + $c->{man_err} = (exists $c->{angle}) ? 35 + $c->{angle} / 2. : 85; + $c->{man_err} = 85 if ($c->{man_err} > 85); # Compare to quaternion used in $m (which provides {ra} {dec} {roll}) which was built - # directly from the 4 components in Backstop - my $q_man = Quat->new($m->{ra}, $m->{dec}, $m->{roll}); - my $q_obc = Quat->new(@c_quat_norm); - my @q_man = @{ $q_man->{q} }; - my $q_diff = $q_man->divide($q_obc); - - if ( abs($q_diff->{ra0} * 3600) > 1.0 - || abs($q_diff->{dec} * 3600) > 1.0 - || abs($q_diff->{roll0} * 3600) > 10.0) - { - push @{ $self->{warn} }, - sprintf( -"Target uplink precision problem for MP_TARGQUAT at $c->{date}\n" - . " Error is yaw, pitch, roll (arcsec) = %.2f %.2f %.2f\n" - . " Use Q1,Q2,Q3,Q4 = %.12f %.12f %.12f %.12f\n", - $q_diff->{ra0} * 3600, - $q_diff->{dec} * 3600, - $q_diff->{roll0} * 3600, - $q_man[0], $q_man[1], $q_man[2], $q_man[3] - ); - } - + # directly from the 4 components in Backstop + my $q_man = Quat->new($m->{ra}, $m->{dec}, $m->{roll}); + my @c_quat_norm = + ($self->{OBC_Q1}, $self->{OBC_Q2}, $self->{OBC_Q3}, $self->{OBC_Q4}); + my $q_obc = Quat->new(@c_quat_norm); + my @q_man = @{ $q_man->{q} }; + my $q_diff = $q_man->divide($q_obc); + + if ( abs($q_diff->{ra0} * 3600) > 1.0 + || abs($q_diff->{dec} * 3600) > 1.0 + || abs($q_diff->{roll0} * 3600) > 10.0) + { + push @{ $self->{warn} }, + sprintf( + "Target uplink precision problem for MP_TARGQUAT at $c->{date}\n" + . " Error is yaw, pitch, roll (arcsec) = %.2f %.2f %.2f\n" + . " Use Q1,Q2,Q3,Q4 = %.12f %.12f %.12f %.12f\n", + $q_diff->{ra0} * 3600, + $q_diff->{dec} * 3600, + $q_diff->{roll0} * 3600, + $q_man[0], $q_man[1], $q_man[2], $q_man[3] + ); } } - push @{ $self->{yellow_warn} }, - sprintf("Did not find match in maneuvers for MP_TARGQUAT at $c->{date}\n") - unless ($found); + + if (not $found) { + push @{ $self->{critical_warn} }, + sprintf( + "Did not find match in maneuvers for MP_TARGQUAT at $c->{date}\n"); + + } } } @@ -2122,10 +2146,10 @@ sub print_report { $o .= sprintf "MP_TARGQUAT at $c->{date} (VCDU count = $c->{vcdu})\n"; $o .= sprintf( " Q1,Q2,Q3,Q4: %.8f %.8f %.8f %.8f\n", - $c->{Q1}, - $c->{Q2}, - $c->{Q3}, - $c->{Q4} + $self->{OBC_Q1}, + $self->{OBC_Q2}, + $self->{OBC_Q3}, + $self->{OBC_Q4} ); if (exists $c->{man_err} and exists $c->{dur} and exists $c->{angle}) { $o .= sprintf( @@ -2995,10 +3019,10 @@ sub proseco_args { obsid => $self->{obsid}, date => $targ_cmd->{stop_date}, att => [ - 0 + $targ_cmd->{q1}, - 0 + $targ_cmd->{q2}, - 0 + $targ_cmd->{q3}, - 0 + $targ_cmd->{q4} + 0 + $self->{OBC_Q1}, + 0 + $self->{OBC_Q2}, + 0 + $self->{OBC_Q3}, + 0 + $self->{OBC_Q4} ], man_angle => 0 + $targ_cmd->{angle}, detector => $si, diff --git a/starcheck/src/starcheck.pl b/starcheck/src/starcheck.pl index f8f7adc1..b9cc4fb4 100755 --- a/starcheck/src/starcheck.pl +++ b/starcheck/src/starcheck.pl @@ -390,11 +390,6 @@ my (%dot_cmd, %dot_time_offset, %dot_tolerance); set_dot_cmd(); -# Go through records and set the time of MP_TARGQUAT commands to -# the time of the subsequent cmd with COMMAND_SW | TLMSID= AOMANUVR - -fix_targquat_time(); - # Now go through records, pull out the interesting things, and assemble # into structures based on obsid. @@ -593,7 +588,7 @@ sub json_obsids { my $obsid_temps; my $json_obsid_temps; $json_obsid_temps = call_python( - "utils.ccd_temp_wrapper", + "calc_ccd_temps.get_ccd_temps", [], { oflsdir => $par{vehicle} ? $STARCHECK : $par{dir}, @@ -735,7 +730,7 @@ sub json_obsids { $out .= "------------ HIGH IR ZONE CHECK -----------------\n\n"; my $ir_report = "${STARCHECK}/high_ir_check.txt"; my $ir_ok = call_python( - "utils.make_ir_check_report", + "check_ir_zone.ir_zone_ok", [], { backstop_file => $backstop, @@ -1100,38 +1095,6 @@ sub make_annotated_file { close $FILE2; } -##*************************************************************************** -sub fix_targquat_time { -##*************************************************************************** - # Go through records and set the time of MP_TARGQUAT commands to - # the time of the subsequent cmd with COMMAND_SW | TLMSID= AOMANUVR - my $manv_time; - my $set = 0; - - for my $i (reverse(0 .. $#cmd)) { - if ($cmd[$i] eq 'COMMAND_SW' and $params[$i] =~ /AOMANUVR/) { - - # print STDERR "First: $cmd[$i], $time[$i], $date[$i] \n"; - $manv_time = $time[$i]; - $set = 1; - } - if ($cmd[$i] eq 'MP_TARGQUAT') { - - # print STDERR "Second: $cmd[$i], $time[$i], $date[$i] \n"; - if ($set eq 1) { - $time[$i] = $manv_time; - - # undef $manv_time; # Make sure that each TARGQUAT gets a unique AOMANUVR time - $set = 0; - } - else { - warning( - "Found MP_TARGQUAT at $date[$i] without corresponding AOMANUVR\n"); - } - } - } -} - ##*************************************************************************** sub set_dot_cmd { ##*************************************************************************** diff --git a/starcheck/utils.py b/starcheck/utils.py index acd797cb..8b1d812b 100644 --- a/starcheck/utils.py +++ b/starcheck/utils.py @@ -30,8 +30,6 @@ import starcheck from starcheck import __version__ as version -from starcheck.calc_ccd_temps import get_ccd_temps -from starcheck.check_ir_zone import ir_zone_ok from starcheck.plot import make_plots_for_obsid ACQS = mica.stats.acq_stats.get_stats() @@ -74,10 +72,6 @@ def time2date(val): return out -def ccd_temp_wrapper(**kwargs): - return get_ccd_temps(**kwargs) - - def plot_cat_wrapper(**kwargs): return make_plots_for_obsid(**kwargs) @@ -99,6 +93,7 @@ def set_kadi_scenario_default(): def get_cheta_source(): + import starcheck.calc_ccd_temps sources = starcheck.calc_ccd_temps.fetch.data_source.sources() if len(sources) == 1 and sources[0] == "cxc": return "cxc" @@ -115,10 +110,6 @@ def get_data_dir(): return sc_data if os.path.exists(sc_data) else "" -def make_ir_check_report(**kwargs): - return ir_zone_ok(**kwargs) - - def get_dither_kadi_state(date): cols = [ "dither", @@ -561,3 +552,19 @@ def vehicle_filter_backstop(backstop_file, outfile): or cmd["type"] == "MP_OBSID"] # Write the filtered commands to the output file write_backstop(filtered_cmds, outfile) + +def target_subfunction_q(q1, q2, q3): + q4_obc = np.sqrt(abs(1.0 - q1**2 - q2**2 - q3**2)) + norm = np.sqrt(q1**2 + q2**2 + q3**2 + q4_obc**2) + q = [q1/norm, q2/norm, q3/norm, q4_obc/norm] + return q + + +def replace_with_obc_quats(bs_cmds): + cmds = bs_cmds.copy() + for cmd in cmds: + if cmd['type'] == 'MP_TARGQUAT': + cmd['params']['q1'], cmd['params']['q2'], cmd['params']['q3'], cmd['params']['q4'] = target_subfunction_q( + cmd['q1'], cmd['q2'], cmd['q3'] + ) + return cmds \ No newline at end of file