#!/usr/bin/env python
#
# Copyright (C) 2010-2015  Jordi Burguet-Castell, Madeline Wade
#
# This program is free software; you can redistribute it and/or modify it
# under the terms of the GNU General Public License as published by the
# Free Software Foundation; either version 2 of the License, or (at your
# option) any later version.
#
# This program is distributed in the hope that it will be useful, but
# WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General
# Public License for more details.
#
# You should have received a copy of the GNU General Public License along
# with this program; if not, write to the Free Software Foundation, Inc.,
# 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.


"""
This pipeline produces h(t) given DARM_ERR and DARM_CTRL or given DELTAL_RESIDUAL and DELTAL_CTRL. It can be run online in real-time or offline on frame files.  It can write h(t) frames to frame files or to a shared memory partition.  

The differential arm length resulting from external sources is

\Delta L_{ext} = d_{err}/(\kappa_c C) + (A_tst * \kappa_tst + A_usum) d_{ctrl}

where C is the sensing function, A_tst is the TST acutuation function, A_usum is the PUM+UIM+TOP actuation, \kappa_c is the time dependent gain of the sensing function, and \kappa_tst in the time-dependent gain of TST actuation.  \Delta L_{ext} is divided by the average arm length (4000 km) to obtain h(t), the external strain in the detectors,

h(t) = \Delta L_{ext} / L .

The time-dependent gains (\kappa's) as well as the value for the coupled cavity pole (f_cc), the time-dependent gain of the PUM actuation (\kappa_pu) and the overall time-depenent gain of the actuation (\kappa_a) are calcuated in this pipeline as well.

This pipeline will most often be run in a format where it picks up after part of the actuation and sensing functions have been applied to the apporiate channels.  In this mode, the input channels are \Delta L_{res} and \Delta L_{ctrl}.  This pipeline then applies further high frequency corrections to each of these channels, applies the appropriate delay to each channel, adds the channels together, and divides by L.

h(t) = (\Delta L_{res} * (1 / \kappa_c) * corrections + (\Delta L_{ctrl, TST} * \kappa_tst + \Delta L_{ctrl, USUM}) * corrections) / L

Note: The \kappa's are complex numbers.  Only the real part of the computed \kappa's are applied as time-dependent gain corrections.

Further documentation explaining the time domain calibration procedure can be found in LIGO DCC #T1400256.

Example: To run this pipeline on some of the earliest locked stretches of the LLO instrument, first query for the locked segments

$ ligolw_segment_query_dqsegdb --query-segments ligolw_segment_query_dqsegdb --segment-url=https://dqsegdb5.phy.syr.edu/ --gps-start-time `lalapps_tconvert November 13 2014` --gps-end-time `lalapps_tconvert now` --include-segments L1:DMT-DC_READOUT --result-name datasegments > L1_segments.xml

The above command will write the locked segments for which the DMT-DC_READOUT flag was on from November 13, 2014 until now to a file called L1_segments.xml. To take a look at the start and stop times 

$ ligolw_print -t segment -c start_time -c end_time L1_segments.xml

Now, use ligo_data_find to find the frame files for these segments,

$ ligo_data_find -o L -t L1_R -s start_time -e end_time -l --url-type file > L1_frames.cache

where you can choose start_time and end_time to be any GPS time within these segments, for example

$ ligo_data_find -o L -t L1_R -s 1100513916 -e 1104641381 -l --url-type file > L1_frames.cache

Now create calibrated data by running this executable with the following command line options with a filters file that is designed to calibration from DARM_ERR and DARM_CTRL

$ gstlal_compute_strain --data-source frames --frame-cache L1_frames.cache --frame-segments-file L1_segments.xml --gps-start-time 1100513916 --gps-end-time 1100514016 --frame-duration=4 --frames-per-file=1 --filters-file FILTERSFILE --frame-segments-name datasegments --ifo L1 --split-actuation-chain --control-sample-rate 16384 --full-calibration

Create calibrated data by running this executable with the following command line options with a filters file that is designed to calibration from DELTAL_CTRL_TST and DELTAL_CTRL_PUM and DELTAL_CTRL_UIM

$ gstlal_compute_strain --data-source frames --frame-cache L1_frames.cache --frame-segments-file L1_segments.xml --gps-start-time 1100513916 --gps-end-time 1100514016 --frame-duration=4 --frames-per-file=1 --filters-file FILTERSFILE --frame-segments-name datasegments --ifo L1 --split-actuation-chain --control-sample-rate 4096 --partial-calibration

The above command requires the filters file, all of which can be found in the calibration SVN. This will write four second frame files containing one four second frame for a just a short, 100 second chunk of time.  You can also write to shared memory by specifying --write-to-shm-partition=PARTITION-NAME.

Type gstlal_compute_strain --help to see the full list of command line options.
"""


import pygtk
pygtk.require("2.0")
import gobject
gobject.threads_init()
import pygst
pygst.require("0.10")
	
import sys
import os
import numpy
import time
import resource

from optparse import OptionParser, Option

# This mess is to make gstreamer stop eating our help messages
if "--help" in sys.argv or "-h" in sys.argv:
	try:
		del sys.argv[sys.argv.index("--help")]
	except ValueError:
		pass
	try:
		del sys.argv[sys.argv.index("-h")]
	except ValueError:
		pass
	import gst
	sys.argv.append("--help")
else:
	import gst

import lal

from gstlal import pipeparts
from gstlal import calibration_parts
from gstlal import simplehandler
from gstlal import datasource

from pylal.xlal.datatypes.ligotimegps import LIGOTimeGPS

from glue.ligolw import ligolw
from glue.ligolw import array
from glue.ligolw import param
from glue.ligolw.utils import segments as ligolw_segments
array.use_in(ligolw.LIGOLWContentHandler)
param.use_in(ligolw.LIGOLWContentHandler)
from glue.ligolw import utils
from glue import segments

def write_graph(demux):
	pipeparts.write_dump_dot(pipeline, "%s.%s" % (options.write_pipeline, "PLAYING"), verbose = True)
	
parser = OptionParser(description = __doc__)

#
# Append program specific options
#

# These options should be used whether the pipeline runs in full calibration mode or partial calibration mode
parser.add_option("--data-source", metavar = "source", help = "Set the data source from [frames|lvshm|white]. Required.")
parser.add_option("--frame-cache", metavar = "filename", help = "Set the name of the LAL cache listing the LIGO .gwf frame files (optional).  This is required iff --data-source=frames")
parser.add_option("--gps-start-time", metavar = "seconds", help = "Set the start time of the segment to analyze in GPS seconds. This is required iff --data-source=frames")
parser.add_option("--gps-end-time", metavar = "seconds", help = "Set the end time of the segment to analyze in GPS seconds. This is required iff --data-source=frames")
parser.add_option("--wings", metavar = "seconds", type = "int", help = "Number of seconds to trim off of the beginning and end of the output. Should only be used if --data-source=frames.")
parser.add_option("--dq-channel-name", metavar = "name", default = "ODC-MASTER_CHANNEL_OUT_DQ", help = "Set the name of the data quality (or state vector) channel. (Default=ODC-MASTER_CHANNEL_OUT_DQ)")
parser.add_option("--ifo", metavar = "name", help = "Name of the IFO to be calibrated.")
parser.add_option("--shared-memory-partition", metavar = "name", help = "Set the name of the shared memory partition to read from.  This is required iff --data-source=lvshm.")
parser.add_option("--frame-segments-file", metavar = "filename", help = "Set the name of the LIGO light-weight XML file from which to load frame segments.  This is required iff --data-source=frames")
parser.add_option("--frame-segments-name", metavar = "name", help = "Set the name of the segments to extract from the segment tables.  This is required iff --frame-segments-file is given")
parser.add_option("--sample-rate", metavar = "Hz", default = 16384, type = "int", help = "Sample rate at which to generate strain data. This should be less than or equal to the sample rate of the error and control signal channels. (Default = 16384 Hz)")
parser.add_option("--control-sample-rate", metavar = "Hz", default = 16384, type = "int", help = "Sample rate of the control signal channels. (Default = 16384 Hz)")
parser.add_option("--odc-sample-rate", metavar = "Hz", default = 16384, type = "int", help = "Sample rate of the ODC state vector channel. (Default = 16384 Hz)")
parser.add_option("--dq-sample-rate", metavar = "Hz", default = 16, type = "int", help = "Sample rate for the outgoing DQ vector GDS-CALIB_STATE_VECTOR. (Default = 16 Hz)")
parser.add_option("--tst-exc-sample-rate", metavar = "Hz", default = 512, type = "int", help = "Sample rate for the control signals being read in. (Default = 512 Hz)")
parser.add_option("--frame-duration", metavar = "seconds", type = "int", default = 4, help = "Set the number of seconds for each frame. (Default = 4)")
parser.add_option("--frames-per-file", metavar = "count", type = "int", default = 1, help = "Set the number of frames per frame file. (Default = 1)")
parser.add_option("--frame-size", metavar = "bytes", type = "int", default = 405338, help = "Approximate size in bytes of frame file images; used when writing to shared memory.  (Default=405338)")
parser.add_option("--compression-scheme", metavar = "scheme", type = "int", default = 256, help = "Set the compression scheme for the framecpp_channelmux element. (Default=256, no compression)")
parser.add_option("--compression-level", metavar = "level", type = "int", default = 0, help = "Set the compression level for the framecpp_channelmux element. (Default=0)")
parser.add_option("--write-to-shm-partition", metavar = "name", help = "Set the name of the shared memory partition to write to. If this is not provided, frames will be written to a file.")
parser.add_option("--buffer-mode", metavar = "number", type = "int", default = 2, help = "Set the buffer mode for the lvshmsink element. (Default=2)")
parser.add_option("--frame-type", metavar = "name", default = "TEST", help = "Set the frame type as input to the frame writing element. (Default=TEST)")
parser.add_option("--output-path", metavar = "name", default = ".", help = "Set the output path for writing frame files. (Default=Current)")
parser.add_option("--no-dq-vector", action = "store_true", help = "Set this if you want to turn off all interpretation and calculation of a data quality vector.")
parser.add_option("--filter-settle-time", metavar = "seconds", type = "int", default = 0, help = "Number of seconds required for the filters to settle in when the data begins to be calibrated. (Default=0)")
parser.add_option("--frequency-domain-filtering", action = "store_true", help = "Set this to perform filtering routines in the frequency domain instead of using direct convolution.")
parser.add_option("--science-quality-bitmask", metavar = "bitmask", type = "int", default = 4, help = "Bitmask used on ODC state vector in order to determine SCIENCE_QUALITY bit information. (Default=4)")
parser.add_option("--science-intent-bitmask", metavar = "bitmask", type = "int", default = 2, help = "Bitmask used on ODC state vector in order to determine SCIENCE_INTENT bit information. (Default=2)")
parser.add_option("--hw-inj-cbc-bitmask", metavar = "bitmask", type = "int", default = 16777216, help = "Bitmask used on ODC state vector in order presence of CBC hardware injection. (Default=16777216)")
parser.add_option("--hw-inj-burst-bitmask", metavar = "bitmask", type = "int", default = 33554432, help = "Bitmask used on ODC state vector in order presence of burst hardware injection. (Default=33554432)")
parser.add_option("--hw-inj-detchar-bitmask", metavar = "bitmask", type = "int", default = 67108864, help = "Bitmask used on ODC state vector in order presence of DetChar hardware injection. (Default=67108864)")
parser.add_option("--hw-inj-stoch-bitmask", metavar = "bitmask", type = "int", default = 8388608, help = "Bitmask used on ODC state vector in order presence of stochastic hardware injection. (Default=8388608)")
parser.add_option("--chan-prefix", metavar = "name", default = "GDS-", help = "Prefix for all output channel names. (Default = GDS)") 
parser.add_option("--chan-suffix", metavar = "name", help = "Suffix for all output channel names.") 

# These are debugging options
parser.add_option("--write-pipeline", metavar = "filename", help = "Write a DOT graph description of the as-built pipeline to this file (optional).  The environment variable GST_DEBUG_DUMP_DOT_DIR must be set for this option to work.")
parser.add_option("-v", "--verbose", action = "store_true", help = "Be verbose (optional).")

# These are options specific to the calibration procedure
parser.add_option("--filters-file", metavar="filename", help = "Name of file containing filters (in npz format)")
parser.add_option("--factors-from-filters-file", action = "store_true", help = "Compute the calibration factors from reference values contained in the filters file instead of from EPICS channels.")
parser.add_option("--no-kappatst", action = "store_true", help = "Set this to turn off the calculation of \kappa_tst.")
parser.add_option("--no-kappapu", action = "store_true", help = "Set this to turn off the calculation of \kappa_pu.")
parser.add_option("--no-kappaa", action = "store_true", help = "Set this to turn off the calculation of \kappa_a.")
parser.add_option("--no-kappac", action = "store_true", help = "Set this to turn off the calculation of \kappa_c.")
parser.add_option("--no-fcc", action = "store_true", help = "Set this to turn off the calculation of f_cc.")
parser.add_option("--factors-averaging-time", metavar = "Sec", type = int, default = 64, help = "Time over which to average the time-varying calibration factors (\kappas), given in seconds. (Default = 64 seconds)")
parser.add_option("--apply-kappaa", action = "store_true", help = "Set this to have the \kappa_a factors multiply the actuation chain.")
parser.add_option("--apply-kappapu", action = "store_true", help = "Set this to have the \kappa_pu factors multiply the actuation chain.")
parser.add_option("--apply-kappatst", action = "store_true", help = "Set this to have the \kappa_tst factors multiply the actuation chain.")
parser.add_option("--apply-kappac", action = "store_true", help = "Set this to have the \kappa_c factors multiply the sensing chain.")
parser.add_option("--compute-factors-sr", metavar = "Hz", type = int, default = 512, help = "Sample rate at which calibration factors are computed. (Default = 512 Hz)")
parser.add_option("--record-factors-sr", metavar = "Hz", type = int, default = 16, help = "Sample rate at which calibration factors are recorded. (Default = 16 Hz)")
parser.add_option("--factors-integration-time", metavar = "sec", type = float, default = 10.0, help = "Integration time for calibration factors computation. (Default = 10 s)")
parser.add_option("--expected-kappaa-real", metavar = "float", type = float, default = 1.0, help = "Expected value for the real part of \kappa_a. (Default = 1.0)")
parser.add_option("--expected-kappapu-real", metavar = "float", type = float, default = 1.0, help = "Expected value for the real part of \kappa_pu. (Default = 1.0)")
parser.add_option("--expected-kappatst-real", metavar = "float", type = float, default = 1.0, help = "Expected value for the real part of \kappa_tst. (Default = 1.0)")
parser.add_option("--expected-kappaa-imag", metavar = "float", type = float, default = 1e-9, help = "Expected value for the imaginary part of \kappa_a. (Default = 0.0)")
parser.add_option("--expected-kappapu-imag", metavar = "float", type = float, default = 1e-9, help = "Expected value for the imaginary part of \kappa_pu. (Default = 0.0)")
parser.add_option("--expected-kappatst-imag", metavar = "float", type = float, default = 1e-9, help = "Expected value for the imaginary part of \kappa_tst. (Default = 0.0)")
parser.add_option("--expected-kappac", metavar = "float", type = float, default = 1.0, help = "Expected value for \kappa_c. (Default = 1.0)")
parser.add_option("--expected-fcc", metavar = "Hz", type = float, default = 330.0, help = "Expected value for the coupled cavity pole. (Default = 330.0 Hz)")
parser.add_option("--kappaa-real-ok-var", metavar = "float", type = float, default = 0.5, help = "Values of the real part of \kappa_a +/- this number will be considered OK. (Default = 0.5)")
parser.add_option("--kappapu-real-ok-var", metavar = "float", type = float, default = 0.5, help = "Values of the real part of \kappa_pu +/- this number will be considered OK. (Default = 0.5)")
parser.add_option("--kappatst-real-ok-var", metavar = "float", type = float, default = 0.5, help = "Values of the real part of \kappa_tst +/- this number will be considered OK. (Default = 0.5)")
parser.add_option("--kappaa-imag-ok-var", metavar = "float", type = float, default = 0.5, help = "Values of the imaginary part of \kappa_a +/- this number will be considered OK. (Default = 0.5)")
parser.add_option("--kappapu-imag-ok-var", metavar = "float", type = float, default = 0.5, help = "Values of the imaginary part of \kappa_pu +/- this number will be considered OK. (Default = 0.5)")
parser.add_option("--kappatst-imag-ok-var", metavar = "float", type = float, default = 0.5, help = "Values of the imaginary part of \kappa_tst +/- this number will be considered OK. (Default = 0.5)")
parser.add_option("--kappac-ok-var", metavar = "float", type = float, default = 0.5, help = "Values of \kappa_c +/- this number will be considered OK. (Default = 0.5)")
parser.add_option("--fcc-ok-var", metavar = "Hz", type = float, default = 30, help = "Values of f_cc +/- this number (in Hz) will be considered OK. (Default = 30 Hz)")
parser.add_option("--exc-channel-name", metavar = "name", default = "CAL-CS_LINE_SUM_DQ", help = "Set the name of the excitation channel.  This is only necessary when the calibration factors computation is turned on, which is the default behavior. (Default = CAL-CS_LINE_SUM_DQ)")
parser.add_option("--tst-exc-channel-name", metavar = "name", default = "SUS-ETMY_L3_CAL_LINE_OUT_DQ", help = "Set the name of the TST excitation channel.  This is only necessary when the \kappa_tst factors computation is turned on, which is the default behavior. (Default = SUS-ETMY_L3_CAL_LINE_OUT_DQ)")
parser.add_option("--pcal-channel-name", metavar = "name", default = "CAL-PCALY_RX_PD_OUT_DQ", help = "Set the name of the PCal channel used for calculating the calibration factors. (Default = CAL-PCALY_RX_PD_OUT_DQ)")
# These are all options related to the reference channels used in the calibraiton factors computation
parser.add_option("--ref-channels-sr", metavar = "Hz", default = 16, help = "Set the sample rate for the reference model channels used in the calibration factors calculation. (Default = 16 Hz)")
parser.add_option("--EP4-real", metavar = "name", default = "CAL-CS_TDEP_DARM_LINE1_REF_A_TST_REAL", help = "Set the name of the channel containing the real part of A_tst at the ESD line used for the \kappa_a and \kappa_pu calculation. (Default = CAL-CS_TDEP_DARM_LINE1_REF_A_TST_REAL)")
parser.add_option("--EP5-real", metavar = "name", default = "CAL-CS_TDEP_DARM_LINE1_REF_A_USUM_REAL", help = "Set the name of the channel containing the real part of A_pu at the ESD line used for the \kappa_a calculation. (Default = CAL-CS_TDEP_DARM_LINE1_REF_A_USUM_REAL)")
parser.add_option("--EP3-real", metavar = "name", default = "CAL-CS_TDEP_DARM_LINE1_REF_A_USUM_INV_REAL", help = "Set the name of the channel containing the real part of 1/A_pu at the ESD line used for the \kappa_pu calculation. (Default = CAL-CS_TDEP_DARM_LINE1_REF_A_USUM_INV_REAL)")
parser.add_option("--EP4-imag", metavar = "name", default = "CAL-CS_TDEP_DARM_LINE1_REF_A_TST_IMAG", help = "Set the name of the channel containing the imaginary part of A_tst at the ESD line used for the \kappa_a and \kappa_pu calculation. (Default = CAL-CS_TDEP_DARM_LINE1_REF_A_TST_IMAG")
parser.add_option("--EP5-imag", metavar = "name", default = "CAL-CS_TDEP_DARM_LINE1_REF_A_USUM_IMAG", help = "Set the name of the channel containing the imaginary part of A_pu at the ESD line used for the \kappa_A calculation. (Default = CAL-CS_TDEP_DARM_LINE1_REF_A_USUM_IMAG")
parser.add_option("--EP3-imag", metavar = "name", default = "CAL-CS_TDEP_DARM_LINE1_REF_A_USUM_INV_IMAG", help = "Set the name of the channel containing the imaginary part of 1/A_pu at the ESD line used for the \kappa_PU calculation. (Default = CAL-CS_TDEP_DARM_LINE1_REF_A_USUM_INV_IMAG")
parser.add_option("--EP2-real", metavar = "name", default = "CAL-CS_TDEP_REF_CLGRATIO_CTRL_REAL", help = "Set the name of the channel containing the real part of the factors used to compute A(f_ctrl). (Default = CAL-CS_TDEP_REF_CLGRATIO_CTRL_REAL)")
parser.add_option("--EP2-imag", metavar = "name", default = "CAL-CS_TDEP_REF_CLGRATIO_CTRL_IMAG", help = "Set the name of the channel containing the imaginary part of the factors used to compute A(f_ctrl). (Default = CAL-CS_TDEP_REF_CLGRATIO_CTRL_IMAG)")
parser.add_option("--EP6-real", metavar = "name", default = "CAL-CS_TDEP_PCALY_LINE2_REF_C_NOCAVPOLE_REAL", help = "Set the name of the channel containing the real part of C_res at the PCal line used for the \kappa_c and f_cc calculation. (Default = CAL-CS_TDEP_PCALY_LINE2_REF_C_NOCAVPOLE_REAL")
parser.add_option("--EP6-imag", metavar = "name", default = "CAL-CS_TDEP_PCALY_LINE2_REF_C_NOCAVPOLE_IMAG", help = "Set the name of the channel containing the imaginary part of C_res at the PCal line used for the \kappa_c and f_cc calculation. (Default = CAL-CS_TDEP_PCALY_LINE2_REF_C_NOCAVPOLE_IMAG")
parser.add_option("--EP7-real", metavar = "name", default = "CAL-CS_TDEP_PCALY_LINE2_REF_D_REAL", help = "Set the name of the channel containing the real part of D at the PCal line used for the \kappa_c and f_cc calculation. (Default = CAL-CS_TDEP_PCALY_LINE2_REF_D_REAL")
parser.add_option("--EP7-imag", metavar = "name", default = "CAL-CS_TDEP_PCALY_LINE2_REF_D_IMAG", help = "Set the name of the channel containing the real part of D at the PCal line used for the \kappa_c and f_cc calculation. (Default = CAL-CS_TDEP_PCALY_LINE2_REF_D_IMAG")
parser.add_option("--EP8-real", metavar = "name", default = "CAL-CS_TDEP_PCALY_LINE2_REF_A_TST_REAL", help = "Set the name of the channel containing the real part of A_tst at the PCal line used for the \kappa_c and f_cc calculation. (Default = CAL-CS_TDEP_PCALY_LINE2_REF_A_TST_REAL")
parser.add_option("--EP8-imag", metavar = "name", default = "CAL-CS_TDEP_PCALY_LINE2_REF_A_TST_IMAG", help = "Set the name of the channel containing the real part of A_tst at the PCal line used for the \kappa_c and f_cc calculation. (Default = CAL-CS_TDEP_PCALY_LINE2_REF_A_TST_IMAG")
parser.add_option("--EP9-real", metavar = "name", default = "CAL-CS_TDEP_PCALY_LINE2_REF_A_USUM_REAL", help = "Set the name of the channel containing the real part of A_pu at the PCal line used for the \kappa_c and f_cc calculation. (Default = CAL-CS_TDEP_PCALY_LINE2_REF_A_USUM_REAL")
parser.add_option("--EP9-imag", metavar = "name", default = "CAL-CS_TDEP_PCALY_LINE2_REF_A_USUM_IMAG", help = "Set the name of the channel containing the real part of A_pu at the PCal line used for the \kappa_c and f_cc calculation. (Default = CAL-CS_TDEP_PCALY_LINE2_REF_A_USUM_IMAG")
parser.add_option("--EP1-real", metavar = "name", default = "CAL-CS_TDEP_REF_INVA_CLGRATIO_TST_REAL", help = "Set the name of the channel containing the real part of the \kappa_tst reference factors. (Default = CAL-CS_TDEP_REF_INVA_CLGRATIO_TST_REAL)")
parser.add_option("--EP1-imag", metavar = "name", default = "CAL-CS_TDEP_REF_INVA_CLGRATIO_TST_IMAG", help = "Set the name of the channel containing the imaginary part of the \kappa_tst reference factors. (Default = CAL-CS_TDEP_REF_INVA_CLGRATIO_TST_IMAG)")
parser.add_option("--split-actuation-chain", action = "store_true", help = "Set this to split the actuation chain into PUM, UIM, TST.")
parser.add_option("--different-control-whitening", action = "store_true", help = "Set when the whitening filters on each section of control chain are different. Only relevant in --split-actuation-chain mode") 

# These options are specific to the full calibration mode
parser.add_option("--full-calibration", action = "store_true", help = "Set this to run the pipeline in full calibration mode.")
parser.add_option("--darm-ctrl-channel-name", metavar = "name", default = "CAL-DARM_CTRL_WHITEN_OUT_DQ", help = "Set the name for the control signal channel. (Default = CAL-DARM_CTRL_WHTIEN_OUT_DQ)")
parser.add_option("--darm-err-channel-name", metavar = "name", default = "CAL-DARM_ERR_WHITEN_OUT_DQ", help = "Set the name of the error signal channel. (Default = CAL-DARM_ERR_WHITEN_OUT_DQ)")
parser.add_option("--control-chain-delay", metavar = "samples", type = int, default = 0, help = "Additional delay from the front end to apply to the control chain. Should be given as number of samples at sample rate specified by --sample-rate (default of 16384 Hz). (Default = 0)")
parser.add_option("--error-chain-delay", metavar = "samples", type = int, default = 0, help = "CALIBRATION MODE: Additional delay from the front end to apply to the error chain. Should be given as number of samples at sample rate specified by --sample-rate (default of 16384 Hz). (Default = 0)")

# These options are specific to the partial calibration mode
parser.add_option("--partial-calibration", action = "store_true", help = "Set this to run the pipeline in partial calibraiton mode.")
parser.add_option("--deltal-ctrl-channel-name", metavar = "name", default = "CAL-DELTAL_CTRL_DQ", help = "Set the name of the partially calibrated control channel. (Default = CAL-DELTAL_CTRL_DQ)")
parser.add_option("--deltal-tst-channel-name", metavar = "name", default = "CAL-DELTAL_CTRL_TST_DQ", help = "Set the name of the partially calibrated control channel for the TST branch of the actuation. (Default = CAL-DELTAL_CTRL_TST_DQ)")
parser.add_option("--deltal-pum-channel-name", metavar = "name", default = "CAL-DELTAL_CTRL_PUM_DQ", help = "Set the name of the partially calibrated control channel for the PUM/UIM branch of the actuation. (Default = CAL-DELTAL_CTRL_PUM_DQ)")
parser.add_option("--deltal-uim-channel-name", metavar = "name", default = "CAL-DELTAL_CTRL_UIM_DQ", help = "Set the name of the partially calibrated control channel for the PUM/UIM branch of the actuation. (Default = CAL-DELTAL_CTRL_UIM_DQ)")
parser.add_option("--deltal-res-channel-name", metavar = "name", default = "CAL-DELTAL_RESIDUAL_DQ", help = "Set the name of the partially calibrated residual channe. (Default = CAL-DELTAL_RESIDUAL_DQ).")

#
# Parse options
#

options, filenames = parser.parse_args()

# Sanity checks for command line options
data_sources = set(("frames", "lvshm"))

if options.data_source not in data_sources:
	raise ValueError("--data-source must be one of %s" % ",".join(data_sources))

if options.data_source == "frames" and options.frame_cache is None:
	raise ValueError("--frame-cache must be specified when using --data-source=frames")

if options.wings is not None and options.data_source != "frames":
	raise ValueError("--wings can only be set when --data-source=frames")

if options.wings is not None and (options.wings % (options.frames_per_file * options.frame_duration)) and (options.wings > (options.frames_per_file * options.frame_duration)):
	raise ValueError("--wings must be an integer multiple of --frames-per-file * --frame-duration when --wings is greater than --frames-per-file * --frame-duration")

if options.wings is not None and (options.wings < (options.frames_per_file * options.frame_duration)) and (((int(options.gps_start_time)+options.wings) % (options.frames_per_file * options.frame_duration) != 0) or ((int(options.gps_end_time)-options.wings) % (options.frames_per_file * options.frame_duration) != 0)):
	raise ValueError("when you want --wings < --frames-per-file * --frame-duration just set your GPS start and end times to be the value of --wings before and after an integer second boundarty of --frames-per-file * --frame-duration")

if options.ifo is None:
	raise ValueError("must specify --ifo")

if options.full_calibration and (options.darm_ctrl_channel_name is None or options.darm_err_channel_name is None):
	raise ValueError("must specify --darm-ctrl-channel-name and --darm-err-channel-name when in full calibration mode") 

if options.partial_calibration and (options.deltal_ctrl_channel_name is None or options.deltal_res_channel_name is None):
	raise ValueError("must specify --deltal-ctrl-channel-name and --deltal-res-channel-name when in partial calibration mode")

if options.partial_calibration and options.split_actuation_chain and (options.deltal_tst_channel_name is None or options.deltal_pum_channel_name is None or options.deltal_uim_channel_name is None or options.deltal_res_channel_name is None):
	raise ValueError("must specify --deltal-tst-channel-name, deltal-pum-channel-name, --deltal-uim-channel-name and deltal-res-channel-name in partial calibration mode when splitting the actuation chain")

if (options.apply_kappatst or options.apply_kappapu) and not options.split_actuation_chain:
	raise ValueError("must specify --split-actuation-chain when applying the \kappa_tst or \kappa_pu correction")

if options.apply_kappaa and (options.apply_kappatst or options.apply_kappapu):
	raise ValueError("cannot apply both \kappa_a and \kappa_tst or \kappa_pu correction as they are not compatible")

if options.frame_segments_file is not None and options.data_source != "frames":
	raise ValueError("can only give --frame-segments-file if --data-source=frames")

if options.frame_segments_name is not None and options.frame_segments_file is None:
	raise ValueError("can only specify --frame-segments-name if --frame-segments-file is given")

if options.data_source == "frames" and (options.gps_start_time is None or options.gps_end_time is None):
	raise ValueError("must specify --gps-start-time and --gps-end-time when --data-source=frames")

if options.full_calibration is None and options.partial_calibration is None:
	raise ValueError("must specify a mode of the pipeline: either --full-calibration or --partial-calibration")

if int(options.record_factors_sr) > int(options.compute_factors_sr):
	raise ValueError("--record-factors-sr must be less than or equal to --compute-factors-sr")

if options.gps_start_time is not None:
	if options.gps_end_time is None:
		raise ValueError("must provide both --gps-start-time and --gps-end-time")
	if options.data_source == "lvshm" or options.data_source == "white":
		raise ValueError("cannot set --gps-start-time or --gps-end-time with --data-source=lvshm or --data-source=white")
	try:
		start = LIGOTimeGPS(options.gps_start_time)
	except ValueError:
		raise ValueError("invalid --gps-start-time %s" % options.gps_start_time)
	try:
		end = LIGOTimeGPS(options.gps_end_time)
	except ValueError:
		raise ValueError("invalid --gps-end-time %s" % options.gps_end_time)
	if start >= end:
		raise ValueError("--gps-start-time must be < --gps-end-time: %s < %s" % (options.gps_start_time, options.gps_end_time))
	# segment from gps start and stop time if given
	seg = segments.segment(start, end)
	# seek event from the gps start and stop time if given
	seekevent = gst.event_new_seek(1., gst.FORMAT_TIME, gst.SEEK_FLAG_FLUSH | gst.SEEK_FLAG_KEY_UNIT, gst.SEEK_TYPE_SET, seg[0].ns(), gst.SEEK_TYPE_SET, seg[1].ns())
elif options.gps_end_time is not None:
	raise ValueError("must provide both --gps-start-time and --gps-end-time")

# Set up instrument and channel name info from command line options
instrument = options.ifo

if options.frame_segments_file is not None:
	# Frame segments from a user defined file
	frame_segments = ligolw_segments.segmenttable_get_by_name(utils.load_filename(options.frame_segments_file, contenthandler = datasource.ContentHandler), options.frame_segments_name).coalesce()
	if seg is not None:
		# clip frame segments to seek segment if it exists (not required, just saves some meory and I/O overhead)
		frame_segments = segments.segmentlistdict((instrument, seglist & segments.segmentlist([seg])) for instrument, seglist in frame_segments.items())
else:
	frame_segments = None

# Set up short-cut names for each of the sample rates used throughout the pipeline and establish caps
sr = options.sample_rate
dqsr = options.dq_sample_rate
odcsr = options.odc_sample_rate
ctrlsr = options.control_sample_rate
caps = "audio/x-raw-float, width=64, rate=%d, channels=1, endianness=1234" % sr # = 8 bytes, a double
ctrl_caps = "audio/x-raw-float, width=64, rate=%d, channels=1, endianness=1234" % ctrlsr 
if options.chan_suffix is not None:
	chan_suffix = options.chan_suffix
else:
	chan_suffix = ""
chan_prefix = options.chan_prefix

td = not options.frequency_domain_filtering

# Make sure we have sufficient resources
# We allocate far more memory than we need, so this is okay

def setrlimit(res, lim):
	hard_lim = resource.getrlimit(res)[1]
	resource.setrlimit(res, (lim if lim is not None else hard_lim, hard_lim))

# set the number of processes and total set size up to hard limit and shrink the per thread stack size (default is 10 MiB)
setrlimit(resource.RLIMIT_NPROC, None)
setrlimit(resource.RLIMIT_AS, None)
setrlimit(resource.RLIMIT_RSS, None)
setrlimit(resource.RLIMIT_STACK, 1024*1024) # 1 Mib per thread

#
# Setup the pipeline
#

pipeline = gst.Pipeline(sys.argv[0])
mainloop = gobject.MainLoop()
handler = simplehandler.Handler(mainloop, pipeline)

# 
# Turn off debugging tools or verboseness
#

pipeparts.mkchecktimestamps = lambda pipeline, src, *args: src # comment this line out to turn on the checktimestamps debugging
if not options.verbose:
	pipeparts.mkprogressreport = lambda pipeline, src, *args: src

#
# Read in data from frames or shared memory
#

if options.data_source == "lvshm":
	src = pipeparts.mklvshmsrc(pipeline, shm_name = options.shared_memory_partition, assumed_duration = 1)
elif options.data_source == "frames":
	src = pipeparts.mklalcachesrc(pipeline, location = options.frame_cache, cache_dsc_regex = instrument)

if options.full_calibration:
	channel_list = [(instrument, options.darm_ctrl_channel_name), (instrument, options.darm_err_channel_name)]
elif options.partial_calibration:
	if options.split_actuation_chain:
		channel_list = [(instrument, options.deltal_tst_channel_name), (instrument, options.deltal_pum_channel_name), (instrument, options.deltal_uim_channel_name), (instrument, options.deltal_res_channel_name)]
	else:
		channel_list = [(instrument, options.deltal_ctrl_channel_name), (instrument, options.deltal_res_channel_name)]
	
if not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappatst or not options.no_kappapu:
	if options.partial_calibration:
		channel_list.append((instrument, options.darm_err_channel_name))
	channel_list.append((instrument, options.pcal_channel_name))

if (not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu) and options.factors_from_filters_file:
	channel_list.append((instrument, options.exc_channel_name))
if (not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu) and not options.factors_from_filters_file:
	channel_list.append((instrument, options.exc_channel_name))
	channel_list.append((instrument, options.EP4_real))
	channel_list.append((instrument, options.EP4_imag))
	if not options.no_kappaa:
		channel_list.append((instrument, options.EP5_real))
		channel_list.append((instrument, options.EP5_imag))
	channel_list.append((instrument, options.EP3_real))
	channel_list.append((instrument, options.EP3_imag))
	channel_list.append((instrument, options.EP2_real))
	channel_list.append((instrument, options.EP2_imag))

if (not options.no_kappac or not options.no_fcc) and not options.factors_from_filters_file:
	channel_list.append((instrument, options.EP6_real))
	channel_list.append((instrument, options.EP6_imag))
	channel_list.append((instrument, options.EP7_real))
	channel_list.append((instrument, options.EP7_imag))
	channel_list.append((instrument, options.EP8_real))
	channel_list.append((instrument, options.EP8_imag))
	channel_list.append((instrument, options.EP9_real))
	channel_list.append((instrument, options.EP9_imag))

if not options.no_kappatst or not options.no_kappaa or not options.no_kappac or not options.no_fcc:
	channel_list.append((instrument, options.tst_exc_channel_name))
	if not options.factors_from_filters_file:
		channel_list.append((instrument, options.EP1_real))
		channel_list.append((instrument, options.EP1_imag))
	
if not options.no_dq_vector:
	channel_list.append((instrument, options.dq_channel_name))
	
# Hook up the relevant channels to the demuxer
if options.data_source == "lvshm":
	demux = pipeparts.mkframecppchanneldemux(pipeline, src, do_file_checksum = True, skip_bad_files = True, channel_list = map("%s:%s".__mod__, channel_list))
elif options.data_source == "frames":
	demux = pipeparts.mkframecppchanneldemux(pipeline, src, do_file_checksum = True, skip_bad_files = False, channel_list = map("%s:%s".__mod__, channel_list))
# Write the pipeline graph after pads have been hooked up to the demuxer
if options.write_pipeline is not None:
	demux.connect("no-more-pads", write_graph)	

# Make sure the code exits when it encounters a non-zero dataValid flag and is running in "offline mode" (i.e. reading from frames)
def exit_on_nonzero_datavalid(pad, ignored):
	datavalid = pad.get_property("datavalid")
	if datavalid != 0:
		print "ERROR: Non-zero dataValid flag encountered"
		sys.exit()

def connect_notify_datavalid_signal(demux, pad):
	name = pad.get_name()
	pad.connect("notify::datavalid", exit_on_nonzero_datavalid)

if options.data_source == "frames":
	demux.connect("pad-added", connect_notify_datavalid_signal)


# Get everything hooked up and blocked off to no more than 1 second buffers
if options.full_calibration:
	ctrl = calibration_parts.hook_up_and_reblock(pipeline, demux, options.darm_ctrl_channel_name, instrument)
	res = calibration_parts.hook_up_and_reblock(pipeline, demux, options.darm_err_channel_name, instrument)
elif options.partial_calibration:
	if options.split_actuation_chain:
		tst = calibration_parts.hook_up_and_reblock(pipeline, demux, options.deltal_tst_channel_name, instrument)
		pum = calibration_parts.hook_up_and_reblock(pipeline, demux, options.deltal_pum_channel_name, instrument)
		uim = calibration_parts.hook_up_and_reblock(pipeline, demux, options.deltal_uim_channel_name, instrument)
	else:
		ctrl = calibration_parts.hook_up_and_reblock(pipeline, demux, options.deltal_ctrl_channel_name, instrument)
	res = calibration_parts.hook_up_and_reblock(pipeline, demux, options.deltal_res_channel_name, instrument)
if not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappatst or not options.no_kappapu:
	pcal = calibration_parts.hook_up_and_reblock(pipeline, demux, options.pcal_channel_name, instrument)
	if options.partial_calibration:
		darm_err = calibration_parts.hook_up_and_reblock(pipeline, demux, options.darm_err_channel_name, instrument)
if (not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu) and options.factors_from_filters_file:
	exc = calibration_parts.hook_up_and_reblock(pipeline, demux, options.exc_channel_name, instrument)
if (not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu) and not options.factors_from_filters_file:
	exc = calibration_parts.hook_up_and_reblock(pipeline, demux, options.exc_channel_name, instrument)
	EP4_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP4_real, instrument)
	EP4_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP4_imag, instrument)
	if not options.no_kappaa:
		EP5_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP5_real, instrument)
		EP5_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP5_imag, instrument)
	EP3_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP3_real, instrument)
	EP3_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP3_imag, instrument)
	EP2_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP2_real, instrument)
	EP2_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP2_imag, instrument)
if (not options.no_kappac or not options.no_fcc) and not options.factors_from_filters_file:
	EP6_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP6_real, instrument)
	EP6_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP6_imag, instrument)
	EP7_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP7_real, instrument)
	EP7_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP7_imag, instrument)
	EP8_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP8_real, instrument)
	EP8_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP8_imag, instrument)
	EP9_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP9_real, instrument)
	EP9_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP9_imag, instrument)
if not options.no_kappatst or not options.no_kappaa or not options.no_kappac or not options.no_fcc:
	tstexc = calibration_parts.hook_up_and_reblock(pipeline, demux, options.tst_exc_channel_name, instrument)
	if not options.factors_from_filters_file:
		EP1_real = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP1_real, instrument)
		EP1_imag = calibration_parts.hook_up_and_reblock(pipeline, demux, options.EP1_imag, instrument)
if not options.no_dq_vector:
	odcstatevector = calibration_parts.hook_up_and_reblock(pipeline, demux, options.dq_channel_name, instrument)
	odcgaptee = pipeparts.mktee(pipeline, odcstatevector)
	odcstatevector = calibration_parts.mkqueue(pipeline, odcgaptee)
	# FIXME: When the ODC is written as unsigned ints, this piece can be removed
	odcstatevector = pipeparts.mkaudioconvert(pipeline, odcstatevector)
	odcstatevector = pipeparts.mkcapsfilter(pipeline, odcstatevector, "audio/x-raw-int, signed=false")

# When reading from disk, clip the incoming data stream(s) to segment list
if options.data_source == "frames" and frame_segments is not None:
	if options.full_calibration or (not options.split_actuation_chain and options.partial_calibration):
		ctrl = pipeparts.mkgate(pipeline, ctrl, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
	elif options.split_actuation_chain and options.partial_calibration:
		tst = pipeparts.mkgate(pipeline, tst, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		pum = pipeparts.mkgate(pipeline, pum, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		uim = pipeparts.mkgate(pipeline, uim, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
	res = pipeparts.mkgate(pipeline, res, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
	if not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappatst or not options.no_kappapu:
		pcal = pipeparts.mkgate(pipeline, pcal, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		if options.partial_calibration:
			darm_err = pipeparts.mkgate(pipeline, darm_err, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
	if (not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu) and options.factors_from_filters_file:
		exc = pipeparts.mkgate(pipeline, exc, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
	if (not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu) and not options.factors_from_filters_file:
		exc = pipeparts.mkgate(pipeline, exc, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP4_real = pipeparts.mkgate(pipeline, EP4_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP4_imag = pipeparts.mkgate(pipeline, EP4_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		if not options.no_kappaa:
			EP5_real = pipeparts.mkgate(pipeline, EP5_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
			EP5_imag = pipeparts.mkgate(pipeline, EP5_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP3_real = pipeparts.mkgate(pipeline, EP3_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP3_imag = pipeparts.mkgate(pipeline, EP3_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP2_real = pipeparts.mkgate(pipeline, EP2_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP2_imag = pipeparts.mkgate(pipeline, EP2_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
	if (not options.no_kappac or not options.no_fcc) and not options.factors_from_filters_file:
		EP6_real = pipeparts.mkgate(pipeline, EP6_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP6_imag = pipeparts.mkgate(pipeline, EP6_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP7_real = pipeparts.mkgate(pipeline, EP7_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP7_imag = pipeparts.mkgate(pipeline, EP7_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP8_real = pipeparts.mkgate(pipeline, EP8_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP8_imag = pipeparts.mkgate(pipeline, EP8_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))	
		EP9_real = pipeparts.mkgate(pipeline, EP9_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		EP9_imag = pipeparts.mkgate(pipeline, EP9_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))	
	if not options.no_kappatst or not options.no_kappaa or not options.no_kappac or not options.no_fcc:
		tstexc = pipeparts.mkgate(pipeline, tstexc, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
		if not options.factors_from_filters_file:
			EP1_real = pipeparts.mkgate(pipeline, EP1_real, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
			EP1_imag = pipeparts.mkgate(pipeline, EP1_imag, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))
	if not options.no_dq_vector:
		odcstatevector = pipeparts.mkgate(pipeline, odcstatevector, threshold = 1, control = pipeparts.mksegmentsrc(pipeline, frame_segments[instrument]))

# Set up the statevector for gating \kappas
if not options.no_dq_vector:
	odcstatevectortee = pipeparts.mktee(pipeline, odcstatevector)
	statevector = pipeparts.mkstatevector(pipeline, calibration_parts.mkqueue(pipeline, odcstatevectortee), required_on = options.science_quality_bitmask)
	statevectortee = pipeparts.mktee(pipeline, statevector)
	odcstatevector = calibration_parts.mkqueue(pipeline, odcstatevectortee)

# Load in the filters file that contains filter coefficients, etc.
filters = numpy.load(options.filters_file)
total_filter_latency = 0

if options.factors_from_filters_file:
	EP1_real = float(filters["EP1_real"])
	EP1_imag = float(filters["EP1_imag"])
	EP2_real = float(filters["EP2_real"])
	EP2_imag = float(filters["EP2_imag"])
	EP3_real = float(filters["EP3_real"])
	EP3_imag = float(filters["EP3_imag"])
	EP4_real = float(filters["EP4_real"])
	EP4_imag = float(filters["EP4_imag"])
	EP5_real = float(filters["EP5_real"])
	EP5_imag = float(filters["EP5_imag"])
	EP6_real = float(filters["EP6_real"])
	EP6_imag = float(filters["EP6_imag"])
	EP7_real = float(filters["EP7_real"])
	EP7_imag = float(filters["EP7_imag"])
	EP8_real = float(filters["EP8_real"])
	EP8_imag = float(filters["EP8_imag"])
	EP9_real = float(filters["EP9_real"])
	EP9_imag = float(filters["EP9_imag"])
if not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappatst or not options.no_kappapu:
	ka_pcal_line_freq = float(filters["ka_pcal_line_freq"])
	ka_pcal_W_real = float(filters["ka_pcal_whitener_re"])
	ka_pcal_W_imag = float(filters["ka_pcal_whitener_im"])
	ka_pcal_corr_real = float(filters["ka_pcal_corr_re"])
	ka_pcal_corr_imag = float(filters["ka_pcal_corr_im"])
if not options.no_kappaa or not options.no_kappac or not options.no_kappatst or not options.no_kappapu:
	ka_esd_line_freq = float(filters["ka_esd_line_freq"])
	ka_esd_W_real = float(filters["ka_esd_whitener_re"])
	ka_esd_W_imag = float(filters["ka_esd_whitener_im"])
if not options.no_kappac or not options.no_fcc:
	kc_pcal_line_freq = float(filters["kc_pcal_line_freq"])
	kc_pcal_W_real = float(filters["kc_pcal_whitener_re"])
	kc_pcal_W_imag = float(filters["kc_pcal_whitener_im"])
	kc_pcal_corr_real = float(filters["kc_pcal_corr_re"])
	kc_pcal_corr_imag = float(filters["kc_pcal_corr_im"])
if not options.no_kappatst or not options.no_kappaa or not options.no_kappac or not options.no_fcc:
	ktst_esd_line_freq = float(filters["ktst_esd_line_freq"])
	ktst_esd_W_real = float(filters["ktst_esd_whitener_re"])
	ktst_esd_W_imag = float(filters["ktst_esd_whitener_im"])
if options.split_actuation_chain:
	if options.partial_calibration:
		tstdewhitensr = int(filters["deltal_tst_dewhiten_sr"])
		pumuimdewhitensr = int(filters["deltal_pumuim_dewhiten_sr"])
		tstdewhitendelay = filters["deltal_tst_dewhiten_delay"]
		pumuimdewhitendelay = filters["deltal_pumuim_dewhiten_delay"]
		tstdewhiten = filters["deltal_tst_dewhiten"]
		pumuimdewhiten = filters["deltal_pumuim_dewhiten"]
		total_filter_latency = total_filter_latency + tstdewhitendelay + pumuimdewhitendelay
	if options.full_calibration:
		tstchainsr = int(filters["actuation_tst_sr"])
		pumuimchainsr = int(filters["actuation_pumuim_sr"])
		tstdelay = filters["actuation_tst_delay"]
		pumuimdelay = filters["actuation_pumuim_delay"]
		tstfilt = filters["actuation_tst"]
		pumuimfilt = filters["actuation_pumuim"]
		ctrldewhitendelay = filters["dewhiten_ctrl_delay"]
		ctrldewhiten = filters["dewhiten_ctrl"]
		ctrldewhitensr = int(filters["dewhiten_ctrl_sr"])
		total_filter_latency = total_filter_latency + pumuimdelay + tstdelay + ctrldewhitendelay
	ctrlcorrsr = int(filters["ctrl_corr_sr"])
	ctrlcorrdelay = filters["ctrl_corr_delay"]
	ctrlcorrfilt = filters["ctrl_corr_filter"]
	total_filter_latency = total_filter_latency + ctrlcorrdelay
elif not options.split_actuation_chain:
	if options.full_calibration:
		ctrlfiltsr = int(filters["actuation_sr"])
		ctrlchaindelay = filters["actuation_delay"]
		ctrlchainfilt = filters["actuation"]
		ctrldewhitendelay = filters["dewhiten_ctrl_delay"]
		ctrldewhiten = filters["dewhiten_ctrl"]
		ctrldewhitensr = int(filters["dewhiten_ctrl_sr"])
		total_filter_latency = total_filter_latency + ctrlchaindelay + ctrldewhitendelay
	if options.partial_calibration:
		ctrldewhitensr = int(filters["deltal_ctrl_dewhiten_sr"])
		ctrldewhitendelay = filters["deltal_ctrl_dewhiten_delay"]
		ctrldewhiten = filters["deltal_ctrl_dewhiten"]
		ctrlfiltsr = int(filters["ctrl_corr_sr"])
		ctrlchaindelay = filters["ctrl_corr_delay"]
		ctrlchainfilt = filters["ctrl_corr_filter"]
		total_filter_latency = total_filter_latency + ctrldewhitendelay + ctrlchaindelay
total_filter_latency_res = 0
if options.full_calibration:
	reschaindelay = filters["inv_sens_delay"]
	reschainfilt = filters["inv_sensing"]
	resdewhitendelay = filters["dewhiten_err_delay"]
	resdewhiten = filters["dewhiten_err"]
elif options.partial_calibration:
	reschaindelay = filters["res_corr_delay"]
	reschainfilt = filters["res_corr_filter"]
	resdewhitendelay = filters["deltal_res_dewhiten_delay"]
	resdewhiten = filters["deltal_res_dewhiten"]
total_filter_latency_res = total_filter_latency_res + reschaindelay + resdewhitendelay

total_filter_latency = max(total_filter_latency, total_filter_latency_res)

#
# TIME-VARYING FACTORS COMPUTATIONS
#

# Set up all of the inputs for the calibration factors
if not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappatst or not options.no_kappapu:
	ref_factors_caps = "audio/x-raw-float, width=64, rate=%d, channels=1, endianness=1234" % options.ref_channels_sr
	compute_calib_factors_caps = "audio/x-raw-float, width=64, rate=%d, channels=1, endianness=1234" % options.compute_factors_sr
	compute_calib_factors_complex_caps = "audio/x-raw-complex, width=128, rate=%d, channels=1, endianness=1234" % options.compute_factors_sr
	pcal = calibration_parts.caps_and_progress(pipeline, pcal, caps, "pcal")
	if options.full_calibration:
		darm_errtee = pipeparts.mktee(pipeline, res)
		darm_err = calibration_parts.mkqueue(pipeline, darm_errtee)
		res = calibration_parts.mkqueue(pipeline, darm_errtee)
	darm_err = calibration_parts.caps_and_progress(pipeline, darm_err, caps, "darm_err")

if (not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu) and options.factors_from_filters_file:
	exc = calibration_parts.caps_and_progress(pipeline, exc, caps, "exc")

if (not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu) and not options.factors_from_filters_file:
	exc = calibration_parts.caps_and_progress(pipeline, exc, caps, "exc")
	EP4_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP4_real, ref_factors_caps, "EP4_real", compute_calib_factors_caps)
	EP4_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP4_imag, ref_factors_caps, "EP4_imag", compute_calib_factors_caps)
	if not options.no_kappaa:
		EP5_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP5_real, ref_factors_caps, "EP5_real", compute_calib_factors_caps)
		EP5_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP5_imag, ref_factors_caps, "EP5_imag", compute_calib_factors_caps)
	EP3_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP3_real, ref_factors_caps, "EP3_real", compute_calib_factors_caps)
	EP3_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP3_imag, ref_factors_caps, "EP3_imag", compute_calib_factors_caps)
	EP2_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP2_real, ref_factors_caps, "EP2_real", compute_calib_factors_caps)
	EP2_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP2_imag, ref_factors_caps, "EP2_imag", compute_calib_factors_caps)
if (not options.no_kappac or not options.no_fcc) and not options.factors_from_filters_file:
	EP6_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP6_real, ref_factors_caps, "EP6_real", compute_calib_factors_caps)
	EP6_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP6_imag, ref_factors_caps, "EP6_imag", compute_calib_factors_caps)
	EP7_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP7_real, ref_factors_caps, "EP7_real", compute_calib_factors_caps)
	EP7_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP7_imag, ref_factors_caps, "EP7_imag", compute_calib_factors_caps)
	EP8_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP8_real, ref_factors_caps, "EP8_real", compute_calib_factors_caps)
	EP8_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP8_imag, ref_factors_caps, "EP8_imag", compute_calib_factors_caps)
	EP9_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP9_real, ref_factors_caps, "EP9_real", compute_calib_factors_caps)
	EP9_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP9_imag, ref_factors_caps, "EP9_imag", compute_calib_factors_caps)
if not options.no_kappatst or not options.no_kappaa or not options.no_kappac or not options.no_fcc:
	tstexccaps = "audio/x-raw-float, rate=%d, channels=1, width=64, endianness=1234" % options.tst_exc_sample_rate
	tstexc = calibration_parts.caps_and_progress(pipeline, tstexc, tstexccaps, "tstexc")
	if not options.factors_from_filters_file:
		EP1_real = calibration_parts.caps_and_progress_and_upsample(pipeline, EP1_real, ref_factors_caps, "EP1_real", compute_calib_factors_caps)
		EP1_imag = calibration_parts.caps_and_progress_and_upsample(pipeline, EP1_imag, ref_factors_caps, "EP1_imag", compute_calib_factors_caps)


# Compute \kappa_a, \kappa_tst,\kappa_c, f_cc, if applicable
if not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappatst or not options.no_kappapu:
	factors_integration_samples = int(options.factors_integration_time) * options.compute_factors_sr
	pcaltee = pipeparts.mktee(pipeline, pcal)
	derrtee = pipeparts.mktee(pipeline, darm_err)
	ka_pcalR_nocorr, ka_pcalI_nocorr = calibration_parts.demodulate(pipeline, calibration_parts.mkqueue(pipeline, pcaltee), sr, ka_pcal_line_freq, caps, compute_calib_factors_caps, factors_integration_samples, td)
	ka_pcalR_nocorr = pipeparts.mkaudioconvert(pipeline, ka_pcalR_nocorr)
	ka_pcalR_nocorr = pipeparts.mkcapsfilter(pipeline, ka_pcalR_nocorr, compute_calib_factors_caps)
	ka_pcalI_nocorr = pipeparts.mkaudioconvert(pipeline, ka_pcalI_nocorr)
	ka_pcalI_nocorr = pipeparts.mkcapsfilter(pipeline, ka_pcalI_nocorr, compute_calib_factors_caps)
	ka_pcalR, ka_pcalI = calibration_parts.filter_at_line(pipeline, calibration_parts.mkqueue(pipeline, ka_pcalR_nocorr), calibration_parts.mkqueue(pipeline, ka_pcalI_nocorr), ka_pcal_corr_real, ka_pcal_corr_imag, compute_calib_factors_caps)
	ka_derrfpR, ka_derrfpI = calibration_parts.demodulate(pipeline, calibration_parts.mkqueue(pipeline, derrtee), sr, ka_pcal_line_freq, caps, compute_calib_factors_caps, factors_integration_samples, td)
	ka_derrfpR = pipeparts.mkaudioconvert(pipeline, ka_derrfpR)
	ka_derrfpR = pipeparts.mkcapsfilter(pipeline, ka_derrfpR, compute_calib_factors_caps)
	ka_derrfpI = pipeparts.mkaudioconvert(pipeline, ka_derrfpI)
	ka_derrfpI = pipeparts.mkcapsfilter(pipeline, ka_derrfpI, compute_calib_factors_caps)
	ka_derrWfpR, ka_derrWfpI = calibration_parts.filter_at_line(pipeline, calibration_parts.mkqueue(pipeline, ka_derrfpR), calibration_parts.mkqueue(pipeline, ka_derrfpI), ka_pcal_W_real, ka_pcal_W_imag, compute_calib_factors_caps)

	pcalfp_derrfpR, pcalfp_derrfpI = calibration_parts.compute_pcalfp_over_derrfp(pipeline, ka_derrWfpR, ka_derrWfpI, ka_pcalR, ka_pcalI, compute_calib_factors_caps)
	pcalfp_derrfpR = pipeparts.mktee(pipeline, pcalfp_derrfpR)
	pcalfp_derrfpI = pipeparts.mktee(pipeline, pcalfp_derrfpI)

if not options.no_kappatst or not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu:
	tstexc = pipeparts.mkresample(pipeline, tstexc, quality = 9)
	tstexc = pipeparts.mkcapsfilter(pipeline, tstexc, caps)
	ktst_tstexcR, ktst_tstexcI = calibration_parts.demodulate(pipeline, tstexc, sr, ktst_esd_line_freq, caps, compute_calib_factors_caps, factors_integration_samples, td)
	ktst_tstexcR = pipeparts.mkaudioconvert(pipeline, ktst_tstexcR)
	ktst_tstexcR = pipeparts.mkcapsfilter(pipeline, ktst_tstexcR, compute_calib_factors_caps)
	ktst_tstexcI = pipeparts.mkaudioconvert(pipeline, ktst_tstexcI)
	ktst_tstexcI = pipeparts.mkcapsfilter(pipeline, ktst_tstexcI, compute_calib_factors_caps)
	ktst_derrftstR, ktst_derrftstI = calibration_parts.demodulate(pipeline, calibration_parts.mkqueue(pipeline, derrtee), sr, ktst_esd_line_freq, caps, compute_calib_factors_caps, factors_integration_samples, td)
	ktst_derrftstR = pipeparts.mkaudioconvert(pipeline, ktst_derrftstR)
	ktst_derrftstR = pipeparts.mkcapsfilter(pipeline, ktst_derrftstR, compute_calib_factors_caps)
	ktst_derrftstI = pipeparts.mkaudioconvert(pipeline, ktst_derrftstI)
	ktst_derrftstI = pipeparts.mkcapsfilter(pipeline, ktst_derrftstI, compute_calib_factors_caps)
	ktst_derrWftstR, ktst_derrWftstI = calibration_parts.filter_at_line(pipeline, calibration_parts.mkqueue(pipeline, ktst_derrftstR), calibration_parts.mkqueue(pipeline, ktst_derrftstI), ktst_esd_W_real, ktst_esd_W_imag, compute_calib_factors_caps)

	if not options.factors_from_filters_file:
		ktstR, ktstI = calibration_parts.compute_kappatst(pipeline, calibration_parts.mkqueue(pipeline, ktst_derrWftstR), calibration_parts.mkqueue(pipeline, ktst_derrWftstI), calibration_parts.mkqueue(pipeline, ktst_tstexcR), calibration_parts.mkqueue(pipeline, ktst_tstexcI), calibration_parts.mkqueue(pipeline, pcalfp_derrfpR), calibration_parts.mkqueue(pipeline, pcalfp_derrfpI), calibration_parts.mkqueue(pipeline, EP1_real), calibration_parts.mkqueue(pipeline, EP1_imag), compute_calib_factors_caps, compute_calib_factors_complex_caps)
	elif options.factors_from_filters_file:
		ktstR, ktstI = calibration_parts.compute_kappatst_from_filters_file(pipeline, calibration_parts.mkqueue(pipeline, ktst_derrWftstR), calibration_parts.mkqueue(pipeline, ktst_derrWftstI), calibration_parts.mkqueue(pipeline, ktst_tstexcR), calibration_parts.mkqueue(pipeline, ktst_tstexcI), calibration_parts.mkqueue(pipeline, pcalfp_derrfpR), calibration_parts.mkqueue(pipeline, pcalfp_derrfpI), EP1_real, EP1_imag, compute_calib_factors_caps, compute_calib_factors_complex_caps)

	ktstRteedq = pipeparts.mktee(pipeline, ktstR)
	ktstIteedq = pipeparts.mktee(pipeline, ktstI)
		 
	if not options.no_dq_vector:
		ktstR = calibration_parts.average_calib_factors(pipeline, calibration_parts.mkqueue(pipeline, ktstRteedq), options.kappatst_real_ok_var, options.expected_kappatst_real, options.factors_averaging_time, compute_calib_factors_caps, 1.0, calibration_parts.mkqueue(pipeline, statevectortee), td)
		ktstI = calibration_parts.average_calib_factors(pipeline, calibration_parts.mkqueue(pipeline, ktstIteedq), options.kappatst_imag_ok_var, options.expected_kappatst_imag, options.factors_averaging_time, compute_calib_factors_caps, 1.0, calibration_parts.mkqueue(pipeline, statevectortee), td)
	ktstRtee = pipeparts.mktee(pipeline, ktstR)
	ktstItee = pipeparts.mktee(pipeline, ktstI)
	if not options.no_kappatst:
		ktstRout = calibration_parts.mkqueue(pipeline, ktstRtee)
		ktstIout = calibration_parts.mkqueue(pipeline, ktstItee)
	
if not options.no_kappaa or not options.no_kappac or not options.no_fcc or not options.no_kappapu:
	ka_excR, ka_excI = calibration_parts.demodulate(pipeline, calibration_parts.mkqueue(pipeline, exc), sr, ka_esd_line_freq,caps, compute_calib_factors_caps, factors_integration_samples, td)
	ka_excR = pipeparts.mkaudioconvert(pipeline, ka_excR)
	ka_excR = pipeparts.mkcapsfilter(pipeline, ka_excR, compute_calib_factors_caps)
	ka_excI = pipeparts.mkaudioconvert(pipeline, ka_excI)
	ka_excI = pipeparts.mkcapsfilter(pipeline, ka_excI, compute_calib_factors_caps)
	ka_derrfxR, ka_derrfxI = calibration_parts.demodulate(pipeline, calibration_parts.mkqueue(pipeline, derrtee), sr, ka_esd_line_freq, caps, compute_calib_factors_caps, factors_integration_samples, td)
	ka_derrfxR = pipeparts.mkaudioconvert(pipeline, ka_derrfxR)
	ka_derrfxR = pipeparts.mkcapsfilter(pipeline, ka_derrfxR, compute_calib_factors_caps)
	ka_derrfxI = pipeparts.mkaudioconvert(pipeline, ka_derrfxI)
	ka_derrfxI = pipeparts.mkcapsfilter(pipeline, ka_derrfxI, compute_calib_factors_caps)
	ka_derrWfxR, ka_derrWfxI = calibration_parts.filter_at_line(pipeline, calibration_parts.mkqueue(pipeline, ka_derrfxR), calibration_parts.mkqueue(pipeline, ka_derrfxI), ka_esd_W_real, ka_esd_W_imag, compute_calib_factors_caps)

	if not options.factors_from_filters_file:
		AfctrlR, AfctrlI = calibration_parts.compute_kappatst(pipeline, ka_derrWfxR, ka_derrWfxI, ka_excR, ka_excI, calibration_parts.mkqueue(pipeline, pcalfp_derrfpR), calibration_parts.mkqueue(pipeline, pcalfp_derrfpI), pipeparts.mkaudioamplify(pipeline, EP2_real, -1.0), pipeparts.mkaudioamplify(pipeline, EP2_imag, -1.0), compute_calib_factors_caps, compute_calib_factors_complex_caps)
	elif options.factors_from_filters_file:
		AfctrlR, AfctrlI = calibration_parts.compute_kappatst_from_filters_file(pipeline, ka_derrWfxR, ka_derrWfxI, ka_excR, ka_excI, calibration_parts.mkqueue(pipeline, pcalfp_derrfpR), calibration_parts.mkqueue(pipeline, pcalfp_derrfpI), -1.0 * EP2_real, -1.0 * EP2_imag, compute_calib_factors_caps, compute_calib_factors_complex_caps)
	AfctrlR = pipeparts.mkaudioconvert(pipeline, AfctrlR)
	AfctrlR = pipeparts.mkcapsfilter(pipeline, AfctrlR, compute_calib_factors_caps)
	AfctrlI = pipeparts.mkaudioconvert(pipeline, AfctrlI)
	AfctrlI = pipeparts.mkcapsfilter(pipeline, AfctrlI, compute_calib_factors_caps)
	if not options.no_kappaa:
		AfctrlR = pipeparts.mktee(pipeline, AfctrlR)
		AfctrlI = pipeparts.mktee(pipeline, AfctrlI)
		if not options.factors_from_filters_file:
			EP4_real = pipeparts.mktee(pipeline, EP4_real)
			EP4_imag = pipeparts.mktee(pipeline, EP4_imag)

	if not options.factors_from_filters_file:
		kpuR, kpuI = calibration_parts.compute_kappapu(pipeline, calibration_parts.mkqueue(pipeline, EP3_real), calibration_parts.mkqueue(pipeline, EP3_imag), calibration_parts.mkqueue(pipeline, AfctrlR), calibration_parts.mkqueue(pipeline, AfctrlI), calibration_parts.mkqueue(pipeline, ktstRtee), calibration_parts.mkqueue(pipeline, ktstItee), calibration_parts.mkqueue(pipeline, EP4_real), calibration_parts.mkqueue(pipeline, EP4_imag), compute_calib_factors_caps, compute_calib_factors_complex_caps)
	elif options.factors_from_filters_file:
		kpuR, kpuI = calibration_parts.compute_kappapu_from_filters_file(pipeline, EP3_real, EP3_imag, calibration_parts.mkqueue(pipeline, AfctrlR), calibration_parts.mkqueue(pipeline, AfctrlI), calibration_parts.mkqueue(pipeline, ktstRtee), calibration_parts.mkqueue(pipeline, ktstItee), EP4_real, EP4_imag, compute_calib_factors_caps, compute_calib_factors_complex_caps)

	kpuRteedq = pipeparts.mktee(pipeline, kpuR)
	kpuIteedq = pipeparts.mktee(pipeline, kpuI)

	if not options.no_dq_vector:
		kpuR = calibration_parts.average_calib_factors(pipeline, calibration_parts.mkqueue(pipeline, kpuRteedq), options.kappapu_real_ok_var, options.expected_kappapu_real, options.factors_averaging_time, compute_calib_factors_caps, 1.0, calibration_parts.mkqueue(pipeline, statevectortee), td)
		kpuI = calibration_parts.average_calib_factors(pipeline, calibration_parts.mkqueue(pipeline, kpuIteedq), options.kappapu_imag_ok_var, options.expected_kappapu_imag, options.factors_averaging_time, compute_calib_factors_caps, 1.0, calibration_parts.mkqueue(pipeline, statevectortee), td)
	kpuRtee = pipeparts.mktee(pipeline, kpuR)
	kpuItee = pipeparts.mktee(pipeline, kpuI)

	if not options.no_kappapu:
		kpuRout = calibration_parts.mkqueue(pipeline, kpuRtee)
		kpuIout = calibration_parts.mkqueue(pipeline, kpuItee)
	
	if not options.no_kappaa:
		if not options.factors_from_filters_file:
			kaR, kaI = calibration_parts.compute_kappaa(pipeline, calibration_parts.mkqueue(pipeline, AfctrlR), calibration_parts.mkqueue(pipeline, AfctrlI), calibration_parts.mkqueue(pipeline, EP4_real), calibration_parts.mkqueue(pipeline, EP4_imag), calibration_parts.mkqueue(pipeline, EP5_real), calibration_parts.mkqueue(pipeline, EP5_imag), compute_calib_factors_caps, compute_calib_factors_complex_caps)
		elif options.factors_from_filters_file:
			kaR, kaI = calibration_parts.compute_kappaa_from_filters_file(pipeline, calibration_parts.mkqueue(pipeline, AfctrlR), calibration_parts.mkqueue(pipeline, AfctrlI), EP4_real, EP4_imag, EP5_real, EP5_imag, compute_calib_factors_caps, compute_calib_factors_complex_caps)

		kaRteedq = pipeparts.mktee(pipeline, kaR)
		kaIteedq = pipeparts.mktee(pipeline, kaI)		

		if not options.no_dq_vector:
			kaR = calibration_parts.average_calib_factors(pipeline, calibration_parts.mkqueue(pipeline, kaRteedq), options.kappaa_real_ok_var, options.expected_kappaa_real, options.factors_averaging_time, compute_calib_factors_caps, 1.0, calibration_parts.mkqueue(pipeline, statevectortee), td)
			kaI = calibration_parts.average_calib_factors(pipeline, calibration_parts.mkqueue(pipeline, kaIteedq), options.kappaa_imag_ok_var, options.expected_kappaa_imag, options.factors_averaging_time, compute_calib_factors_caps, 1.0, calibration_parts.mkqueue(pipeline, statevectortee), td)
		kaRtee = pipeparts.mktee(pipeline, kaR)
		kaItee = pipeparts.mktee(pipeline, kaI)
		kaRout = calibration_parts.mkqueue(pipeline, kaRtee)
		kaIout = calibration_parts.mkqueue(pipeline, kaItee)

	if not options.no_kappac or not options.no_fcc:
		kc_pcalR_nocorr, kc_pcalI_nocorr = calibration_parts.demodulate(pipeline, calibration_parts.mkqueue(pipeline, pcaltee), sr, kc_pcal_line_freq, caps, compute_calib_factors_caps, factors_integration_samples, td)
		kc_pcalR_nocorr = pipeparts.mkaudioconvert(pipeline, kc_pcalR_nocorr)
		kc_pcalR_nocorr = pipeparts.mkcapsfilter(pipeline, kc_pcalR_nocorr, compute_calib_factors_caps)
		kc_pcalI_nocorr = pipeparts.mkaudioconvert(pipeline, kc_pcalI_nocorr)
		kc_pcalI_nocorr = pipeparts.mkcapsfilter(pipeline, kc_pcalI_nocorr, compute_calib_factors_caps)
		kc_pcalR, kc_pcalI = calibration_parts.filter_at_line(pipeline, calibration_parts.mkqueue(pipeline, kc_pcalR_nocorr), calibration_parts.mkqueue(pipeline, kc_pcalI_nocorr), kc_pcal_corr_real, kc_pcal_corr_imag, compute_calib_factors_caps)
		kc_derrR, kc_derrI = calibration_parts.demodulate(pipeline, calibration_parts.mkqueue(pipeline, derrtee), sr, kc_pcal_line_freq, caps, compute_calib_factors_caps, factors_integration_samples, td)
		kc_derrR = pipeparts.mkaudioconvert(pipeline, kc_derrR)
		kc_derrR = pipeparts.mkcapsfilter(pipeline, kc_derrR, compute_calib_factors_caps)
		kc_derrI = pipeparts.mkaudioconvert(pipeline, kc_derrI)
		kc_derrI = pipeparts.mkcapsfilter(pipeline, kc_derrI, compute_calib_factors_caps)
		kc_derrWR, kc_derrWI = calibration_parts.filter_at_line(pipeline, calibration_parts.mkqueue(pipeline, kc_derrR), calibration_parts.mkqueue(pipeline, kc_derrI), kc_pcal_W_real, kc_pcal_W_imag, compute_calib_factors_caps)
		kc_pcal_derrR, kc_pcal_derrI = calibration_parts.compute_pcalfp_over_derrfp(pipeline, kc_derrWR, kc_derrWI, kc_pcalR, kc_pcalI, compute_calib_factors_caps)
		
		if not options.factors_from_filters_file:
			SR, SI = calibration_parts.compute_S(pipeline, calibration_parts.mkqueue(pipeline, EP6_real), calibration_parts.mkqueue(pipeline, EP6_imag), calibration_parts.mkqueue(pipeline, kc_pcal_derrR), calibration_parts.mkqueue(pipeline, kc_pcal_derrI), calibration_parts.mkqueue(pipeline, ktstRtee), calibration_parts.mkqueue(pipeline, ktstItee), calibration_parts.mkqueue(pipeline, kpuRtee), calibration_parts.mkqueue(pipeline, kpuItee), calibration_parts.mkqueue(pipeline, EP7_real), calibration_parts.mkqueue(pipeline, EP7_imag), calibration_parts.mkqueue(pipeline, EP8_real), calibration_parts.mkqueue(pipeline, EP8_imag), calibration_parts.mkqueue(pipeline, EP9_real), calibration_parts.mkqueue(pipeline, EP9_imag), compute_calib_factors_caps, compute_calib_factors_complex_caps)
		elif options.factors_from_filters_file:
			SR, SI = calibration_parts.compute_S_from_filters_file(pipeline, EP6_real, EP6_imag, calibration_parts.mkqueue(pipeline, kc_pcal_derrR), calibration_parts.mkqueue(pipeline, kc_pcal_derrI), calibration_parts.mkqueue(pipeline, ktstRtee), calibration_parts.mkqueue(pipeline, ktstItee), calibration_parts.mkqueue(pipeline, kpuRtee), calibration_parts.mkqueue(pipeline, kpuItee), EP7_real, EP7_imag, EP8_real, EP8_imag, EP9_real, EP9_imag, compute_calib_factors_caps, compute_calib_factors_complex_caps)
		SR = pipeparts.mkaudioconvert(pipeline, SR)
		SR = pipeparts.mkcapsfilter(pipeline, SR, compute_calib_factors_caps)
		SRtee = pipeparts.mktee(pipeline, SR)
		SI = pipeparts.mkaudioconvert(pipeline, SI)
		SI = pipeparts.mkcapsfilter(pipeline, SI, compute_calib_factors_caps)
		SItee = pipeparts.mktee(pipeline, SI)
		if not options.no_kappac:
			kc = calibration_parts.compute_kappac(pipeline, calibration_parts.mkqueue(pipeline, SRtee), calibration_parts.mkqueue(pipeline, SItee), compute_calib_factors_caps)
			
			kcteedq = pipeparts.mktee(pipeline, kc)

			if not options.no_dq_vector:
				kc = calibration_parts.average_calib_factors(pipeline, calibration_parts.mkqueue(pipeline, kcteedq), options.kappac_ok_var, options.expected_kappac, options.factors_averaging_time, compute_calib_factors_caps, 1.0, calibration_parts.mkqueue(pipeline, statevectortee), td)
			kctee = pipeparts.mktee(pipeline, kc)
			kcout = calibration_parts.mkqueue(pipeline, kctee)
		if not options.no_fcc:
			fcc = calibration_parts.compute_fcc(pipeline, calibration_parts.mkqueue(pipeline, SRtee), calibration_parts.mkqueue(pipeline, SItee), kc_pcal_line_freq, compute_calib_factors_caps)
			fccteedq = pipeparts.mktee(pipeline, fcc)
			if not options.no_dq_vector:
				fcc = calibration_parts.average_calib_factors(pipeline, calibration_parts.mkqueue(pipeline, fccteedq), options.fcc_ok_var, options.expected_fcc, options.factors_averaging_time, compute_calib_factors_caps, 1.0, calibration_parts.mkqueue(pipeline, statevectortee), td)
			fcctee = pipeparts.mktee(pipeline, fcc)
			fccout = calibration_parts.mkqueue(pipeline, fcctee)
#
# CONTROL BRANCH
#

# The reverse of the filters will be used in all filtering below due to the definition of the filtering procedure employed by lal_firbank
if options.split_actuation_chain:
	ctrlcorrcaps = "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % ctrlcorrsr
	if options.partial_calibration:
		tst = calibration_parts.caps_and_progress(pipeline, tst, ctrl_caps, "tst")
		tst = pipeparts.mkfirbank(pipeline, tst, fir_matrix = [[0,1]], time_domain = td)
		if options.different_control_whitening:
			tst = calibration_parts.resample(pipeline, tst, "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % tstdewhitensr)
			tst = pipeparts.mkfirbank(pipeline, tst, latency = int(tstdewhitendelay), fir_matrix = [tstdewhiten[::-1]], time_domain = td)
		tst = calibration_parts.resample(pipeline, tst, ctrlcorrcaps)
	
		pum = calibration_parts.caps_and_progress(pipeline, pum, ctrl_caps, "pum")
		uim = calibration_parts.caps_and_progress(pipeline, uim, ctrl_caps, "uim")
		pumuim = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, pum, uim), ctrl_caps)
		pumuim = pipeparts.mkaudioconvert(pipeline, pumuim)
		pumuim = pipeparts.mkcapsfilter(pipeline, pumuim, ctrl_caps) 
		pumuim = pipeparts.mkfirbank(pipeline, pumuim, fir_matrix = [[0,1]], time_domain = td)
		if options.different_control_whitening:
			pumuim = calibration_parts.resample(pipeline, pumuim, "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % pumuimdewhitensr) 
			pumuim = pipeparts.mkfirbank(pipeline, pumuim, latency = int(pumuimdewhitendelay), fir_matrix = [pumuimdewhiten[::-1]], time_domain = td)
		pumuim = calibration_parts.resample(pipeline, pumuim, ctrlcorrcaps)

	if options.full_calibration:
		ctrl = calibration_parts.caps_and_progress(pipeline, ctrl, caps, "ctrl")
		ctrl = pipeparts.mkfirbank(pipeline, ctrl, fir_matrix = [[0,1]], time_domain = td)
		ctrl = calibration_parts.resample(pipeline, ctrl, "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % ctrldewhitensr)
		ctrl = pipeparts.mkfirbank(pipeline, ctrl, latency = int(ctrldewhitendelay), fir_matrix = [ctrldewhiten[::-1]], time_domain = td)
		ctrltee = pipeparts.mktee(pipeline, ctrl)

		tst = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, ctrltee), "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" %  tstchainsr)
		tst = pipeparts.mkfirbank(pipeline, tst, latency = int(tstdelay), fir_matrix = [tstfilt[::-1]], time_domain = td)
		tst = calibration_parts.resample(pipeline, tst, ctrlcorrcaps)

		pumuim = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, ctrltee), "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % pumuimchainsr)
		pumuim = pipeparts.mkfirbank(pipeline, pumuim, latency = int(pumuimdelay), fir_matrix = [pumuimfilt[::-1]], time_domain = td)
		pumuim = calibration_parts.resample(pipeline, pumuim, ctrlcorrcaps)

	if options.apply_kappatst:
		# Only apply the real part of \kappa_tst as a correction to A_tst
		ktst_for_tst = calibration_parts.mkqueue(pipeline, ktstRtee)
		ktst_for_tst = calibration_parts.resample(pipeline, ktst_for_tst, ctrlcorrcaps)
		ktst_for_tst = pipeparts.mkgeneric(pipeline, ktst_for_tst, "lal_check_calib_factors", min=options.expected_kappatst_real-options.kappatst_real_ok_var, max=options.expected_kappatst_real+options.kappatst_real_ok_var, default = options.expected_kappatst_real)
		tst = calibration_parts.mkmultiplier(pipeline, calibration_parts.list_srcs(pipeline, ktst_for_tst, tst), ctrlcorrcaps)
		tst = pipeparts.mkaudioconvert(pipeline, tst)
		tst = pipeparts.mkcapsfilter(pipeline, tst, ctrlcorrcaps)
	if options.apply_kappapu:
		# Only apply the real part of \kappa_pu as a correction to A_pu
		kpu_for_pu = calibration_parts.mkqueue(pipeline, kpuRtee)
		kpu_for_pu = calibration_parts.resample(pipeline, kpu_for_pu, ctrlcorrcaps)
		kpu_for_pu = pipeparts.mkgeneric(pipeline, kpu_for_pu, "lal_check_calib_factors", min=options.expected_kappapu_real-options.kappapu_real_ok_var, max=options.expected_kappapu_real+options.kappapu_real_ok_var, default = options.expected_kappapu_real)
		pumuim = calibration_parts.mkmultiplier(pipeline, calibration_parts.list_srcs(pipeline, kpu_for_pu, pumuim), ctrlcorrcaps)
		pumuim = pipeparts.mkaudioconvert(pipeline, pumuim)
		pumuim = pipeparts.mkcapsfilter(pipeline, pumuim, ctrlcorrcaps)

	ctrl = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, tst, pumuim), ctrlcorrcaps)
	ctrl = pipeparts.mkaudioconvert(pipeline, ctrl)
	ctrl = pipeparts.mkcapsfilter(pipeline, ctrl, ctrlcorrcaps)
	if options.apply_kappaa:
		ka_modify_ctrl = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, kaRtee), caps)
		ka_modify_ctrl = pipeparts.mkgeneric(pipeline, ka_modify_ctrl, "lal_check_calib_factors", min=options.expected_kappaa_real-options.kappaa_real_ok_var, max=options.expected_kappaa_real+options.kappaa_real_ok_var, default = options.expected_kappaa_real)
		ctrl = calibration_parts.mkmultiplier(pipeline, calibration_parts.list_srcs(pipeline, ctrl, ka_modify_ctrl), ctrlcorrcaps)
		ctrl = pipeparts.mkaudioconvert(pipeline, ctrl)
		ctrl = pipeparts.mkcapsfilter(pipeline, ctrl, ctrlcorrcaps)
	if not options.different_control_whitening and options.partial_calibration:
		ctrl = calibration_parts.resample(pipeline, ctrl, "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % pumuimdewhitensr) 
		ctrl = pipeparts.mkfirbank(pipeline, ctrl, latency = int(pumuimdewhitendelay), fir_matrix = [pumuimdewhiten[::-1]], time_domain = td)
		ctrl = calibration_parts.resample(pipeline, ctrl, ctrlcorrcaps)
	if options.partial_calibration:
		ctrl = pipeparts.mkfirbank(pipeline, ctrl, latency = int(ctrlcorrdelay), fir_matrix = [ctrlcorrfilt[::-1]], time_domain = td)
elif not options.split_actuation_chain:
	ctrl = calibration_parts.caps_and_progress(pipeline, ctrl, ctrl_caps, "ctrl")
	ctrl = pipeparts.mkfirbank(pipeline, ctrl, fir_matrix = [[0,1]], time_domain = td)
	ctrl = calibration_parts.resample(pipeline, ctrl, "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % ctrlfiltsr)
	ctrl = pipeparts.mkfirbank(pipeline, ctrl, latency = int(ctrlchaindelay), fir_matrix = [ctrlchainfilt[::-1]], time_domain = td)
	ctrl = calibration_parts.resample(pipeline, ctrl, "audio/x-raw-float, width=64, endianness=1234, channels=1, rate=%d" % ctrldewhitensr)
	ctrl = pipeparts.mkfirbank(pipeline, ctrl, latency = int(ctrldewhitendelay), fir_matrix = [ctrldewhiten[::-1]], time_domain = td)
ctrl = calibration_parts.resample(pipeline, ctrl, caps)

#
# RESIDUAL BRANCH
#

res = calibration_parts.caps_and_progress(pipeline, res, caps, "res")
res = pipeparts.mkfirbank(pipeline, res, fir_matrix = [[0,1]], time_domain = td) # FIXME: Not sure why this is needed... test more without it or debug the underlying issue.

# Apply factors to actuation and sensing chains, if applicable
if options.apply_kappac:
	kc_modify_res = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, kctee), "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % sr)
	kc_modify_res = pipeparts.mkgeneric(pipeline, kc_modify_res, "lal_check_calib_factors", min=options.expected_kappac-options.kappac_ok_var, max=options.expected_kappac+options.kappac_ok_var, default = options.expected_kappac)
	res = calibration_parts.mkmultiplier(pipeline, calibration_parts.list_srcs(pipeline, res, pipeparts.mkpow(pipeline, kc_modify_res, exponent = -1.0)), "audio/x-raw-float, width=64, rate=%d, channels=1, endianness=1234" % sr)
	res = pipeparts.mkaudioconvert(pipeline, res)
	res = pipeparts.mkcapsfilter(pipeline, res, "audio/x-raw-float, width=64, rate=%d, channels=1, endianness=1234" % sr)

# The reverse of the filters will be used in all filtering below due to the definition of the filtering procedure employed by lal_firbank
res = pipeparts.mkfirbank(pipeline, res, latency = int(reschaindelay), fir_matrix = [reschainfilt[::-1]], time_domain = td)
res = pipeparts.mkfirbank(pipeline, res, latency = int(resdewhitendelay), fir_matrix = [resdewhiten[::-1]], time_domain = td)
		
#
# CONTROL + RESIDUAL = H(T)
#

# Add control and residual chains and divide by L to make h(t)
strain = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, res, ctrl), caps)
# Divide by L in a way that is compatitble with old and new filters files, since old filter files don't recored "arm length"
try:
	strain = pipeparts.mkaudioamplify(pipeline, strain, 1.0/float(filters["arm_length"]))
except KeyError:
	strain = pipeparts.mkaudioamplify(pipeline, strain, 1.0/3995.1)
strain = calibration_parts.mkaudiorate(pipeline, strain)
strain = calibration_parts.mkreblock(pipeline, strain) 
strain = pipeparts.mkprogressreport(pipeline, strain, "progress_hoft_%s" % instrument)
	
# Put the units back to strain before writing to frames
straintee = pipeparts.mktee(pipeline, strain)
straintagstr = "units=strain,channel-name=%sCALIB_STRAIN%s,instrument=%s" % (chan_prefix, chan_suffix, instrument)
strain = pipeparts.mktaginject(pipeline, calibration_parts.mkqueue(pipeline, straintee), straintagstr)

#
# GDS-CALIB_STATE_VECTOR BRANCH
#

if not options.no_dq_vector:
	odcstatevector = calibration_parts.mkaudiorate(pipeline, odcstatevector)
	odcstatevector = calibration_parts.mkreblock(pipeline, odcstatevector) 
	odctagstr = "channel-name=%s:%s, instrument=%s" % (instrument, options.dq_channel_name, instrument)
	odcstatevector = pipeparts.mktaginject(pipeline, odcstatevector, odctagstr)
	odcstatevector = pipeparts.mkprogressreport(pipeline, odcstatevector, "progress_odc_%s" % instrument)
	odcstatevectortee = pipeparts.mktee(pipeline, odcstatevector)

	#
	# GAP BIT BRANCH
	#	

	nogap = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, odcstatevectortee), control = calibration_parts.mkqueue(pipeline, odcgaptee), threshold = 0)
	nogap = pipeparts.mkbitvectorgen(pipeline, nogap, bit_vector = 16384, nongap_is_control = True)
	nogap = pipeparts.mkcapsfilter(pipeline, nogap, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, rate=%d, endianness=1234" % odcsr)
	nogap = calibration_parts.mkaudiorate(pipeline, nogap)
	nogap = calibration_parts.mkreblock(pipeline, nogap)
	nogap = pipeparts.mkgeneric(pipeline, nogap, "lal_logical_undersampler", required_on = 16384, status_out = 16384)
	nogap = pipeparts.mkcapsfilter(pipeline, nogap, "audio/x-raw-int, rate=%d, width=32, depth=32, signed=false, channels=1, endianness=1234" % dqsr)

	# 
	# SCIENCE-INTENT BIT BRANCH
	#

	scienceintent = calibration_parts.mkqueue(pipeline, odcstatevectortee)
	scienceintent = pipeparts.mkgeneric(pipeline, scienceintent, "lal_logical_undersampler", required_on = options.science_intent_bitmask, status_out = 2)
	scienceintent = pipeparts.mkcapsfilter(pipeline, scienceintent, "audio/x-raw-int, width=32, depth=32, endianness=1234, channels=1, signed=false, rate=%d" % dqsr)
	scienceintent = calibration_parts.mkaudiorate(pipeline, scienceintent)
	scienceintent = calibration_parts.mkreblock(pipeline, scienceintent)
	scienceintenttee = pipeparts.mktee(pipeline, scienceintent)	

	#
	# SCIENCE-QUALITY BIT BRANCH
	#

	sciencequality = pipeparts.mkgeneric(pipeline, calibration_parts.mkqueue(pipeline, odcstatevectortee), "lal_logical_undersampler", required_on = options.science_quality_bitmask, status_out = 4)
	sciencequality = pipeparts.mkcapsfilter(pipeline, sciencequality, "audio/x-raw-int, width=32, depth=32, channels=1, endianness=1234, signed=false, rate=%d" % dqsr)
	sciencequality = calibration_parts.mkaudiorate(pipeline, sciencequality)
	sciencequality = calibration_parts.mkreblock(pipeline, sciencequality)
	sciencequalitytee = pipeparts.mktee(pipeline, sciencequality)
	
	#
	# H(t)-PRODUCED BIT BRANCH
	#

	htproduced = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, straintee), control = calibration_parts.mkqueue(pipeline, straintee), threshold = 0)
	htproduced = pipeparts.mkbitvectorgen(pipeline, htproduced, bit_vector = 8, nongap_is_control = True)
	htproduced = pipeparts.mkcapsfilter(pipeline, htproduced, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, rate=%d, endianness=1234" % sr)
	htproduced = pipeparts.mkgeneric(pipeline, htproduced, "lal_logical_undersampler", required_on = 8, status_out = 8)
	htproduced = pipeparts.mkcapsfilter(pipeline, htproduced, "audio/x-raw-int, rate=%d, width=32, depth=32, signed=false, channels=1, endianness=1234" % dqsr)
	htproduced = calibration_parts.mkaudiorate(pipeline, htproduced)
	htproduced = calibration_parts.mkreblock(pipeline, htproduced)

	#
	# FILTERS-OK BIT BRANCH
	#
	
	# Set the FILTERS-OK bit based on science-quality transitions
	filtersok = pipeparts.mkbitvectorgen(pipeline, calibration_parts.mkqueue(pipeline, sciencequalitytee), bit_vector=16, threshold=4)
	filtersok = pipeparts.mkcapsfilter(pipeline, filtersok, "audio/x-raw-int, width=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)
	filtersok = calibration_parts.mkaudiorate(pipeline, filtersok)
	filtersok = calibration_parts.mkreblock(pipeline, filtersok)
	filtersok = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, filtersok), control = calibration_parts.mkqueue(pipeline, sciencequalitytee), threshold = 4, attack_length = -int(options.filter_settle_time) * dqsr)

	#
	# KAPPAA-OK BIT BRANCH
	#
	if not options.no_kappaa:
		kaIdq = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, kaIteedq), "audio/x-raw-float, rate=%d, width=64, channels=1, endianness=1234" % dqsr)
		kaIdq = pipeparts.mkgeneric(pipeline, kaIdq, "lal_add_constant", constant = -options.expected_kappaa_imag)
		kaRdq = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, kaRteedq), "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % dqsr)
		kaRdq = pipeparts.mkgeneric(pipeline, kaRdq, "lal_add_constant", constant = -options.expected_kappaa_real)
		
		kaok = pipeparts.mkbitvectorgen(pipeline, kaRdq, threshold = options.kappaa_real_ok_var, invert_control = True, bit_vector = 512)
		kaok = pipeparts.mkcapsfilter(pipeline, kaok, "audio/x-raw-int, depth=32, width=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)
		kaok = calibration_parts.mkaudiorate(pipeline, kaok)
		kaok = calibration_parts.mkreblock(pipeline, kaok)
		kaok = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, kaok), threshold = options.kappaa_imag_ok_var, invert_control = True, control = calibration_parts.mkqueue(pipeline, kaIdq))
		kaok = pipeparts.mkcapsfilter(pipeline, kaok, "audio/x-raw-int, depth=32, width=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)

	#
	# KAPPAPU-OK BIT BRANCH
	#
	if not options.no_kappapu:
		kpuIdq = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, kpuIteedq), "audio/x-raw-float, rate=%d, width=64, channels=1, endianness=1234" % dqsr)
		kpuIdq = pipeparts.mkgeneric(pipeline, kpuIdq, "lal_add_constant", constant = -options.expected_kappapu_imag)
		kpuRdq = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, kpuRteedq), "audio/x-raw-float, width=64, channels=1, endianness=1234, rate=%d" % dqsr)
		kpuRdq = pipeparts.mkgeneric(pipeline, kpuRdq, "lal_add_constant", constant = -options.expected_kappapu_real)
		
		kpuok = pipeparts.mkbitvectorgen(pipeline, kpuRdq, threshold = options.kappapu_real_ok_var, invert_control = True, bit_vector = 1024)
		kpuok = pipeparts.mkcapsfilter(pipeline, kpuok, "audio/x-raw-int, depth=32, width=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)
		kpuok = calibration_parts.mkaudiorate(pipeline, kpuok)
		kpuok = calibration_parts.mkreblock(pipeline, kpuok)
		kpuok = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, kpuok), threshold = options.kappapu_imag_ok_var, invert_control = True, control = calibration_parts.mkqueue(pipeline, kpuIdq))
		kpuok = pipeparts.mkcapsfilter(pipeline, kpuok, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)

	#
	# KAPPATST-OK BIT BRANCH
	#
	if not options.no_kappatst:
		ktstIdq = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, ktstIteedq), "audio/x-raw-float, rate=%d, width=64, channels=1, endianness=1234" % dqsr)
		ktstIdq = pipeparts.mkgeneric(pipeline, ktstIdq, "lal_add_constant", constant = -options.expected_kappatst_imag)
		ktstRdq = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, ktstRteedq), "audio/x-raw-float, rate=%d, width=64, channels=1, endianness=1234" % dqsr)
		ktstRdq = pipeparts.mkgeneric(pipeline, ktstRdq, "lal_add_constant", constant = -options.expected_kappatst_real)
		
		ktstok = pipeparts.mkbitvectorgen(pipeline, ktstRdq, threshold = options.kappatst_real_ok_var, invert_control = True, bit_vector = 2048)
		ktstok = pipeparts.mkcapsfilter(pipeline, ktstok, "audio/x-raw-int, depth=32, width=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)
		ktstok = calibration_parts.mkaudiorate(pipeline, ktstok)
		ktstok = calibration_parts.mkreblock(pipeline, ktstok)
		ktstok = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, ktstok), threshold = options.kappatst_imag_ok_var, invert_control = True, control = calibration_parts.mkqueue(pipeline, ktstIdq)) 
		ktstok = pipeparts.mkcapsfilter(pipeline, ktstok, "audio/x-raw-int, depth=32, width=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)

	#
	# KAPPAC-OK BIT BRANCH
	#
	if not options.no_kappac:
		kcdq = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, kcteedq), "audio/x-raw-float, rate=%d, width=64, endianness=1234, channels=1" % dqsr)
		kcdq = pipeparts.mkgeneric(pipeline, kcdq, "lal_add_constant", constant = -options.expected_kappac)
		
		kcok = pipeparts.mkbitvectorgen(pipeline, kcdq, threshold = options.kappac_ok_var, invert_control = True, bit_vector = 4096)
		kcok = pipeparts.mkcapsfilter(pipeline, kcok, "audio/x-raw-int, depth=32, width=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)
		kcok = calibration_parts.mkaudiorate(pipeline, kcok)
		kcok = calibration_parts.mkreblock(pipeline, kcok)

	#
	# FCC-OK BIT BRANCH
	#
	if not options.no_fcc:
		fccdq = calibration_parts.resample(pipeline, calibration_parts.mkqueue(pipeline, fccteedq), "audio/x-raw-float, rate=%d, width=64, channels=1, endianness=1234" % dqsr)
		fccdq = pipeparts.mkgeneric(pipeline, fccdq, "lal_add_constant", constant = -options.expected_fcc)
		
		fccok = pipeparts.mkbitvectorgen(pipeline, fccdq, threshold = options.fcc_ok_var, invert_control = True, bit_vector = 8192)
		fccok = pipeparts.mkcapsfilter(pipeline, fccok, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)
		fccok = calibration_parts.mkaudiorate(pipeline, fccok)
		fccok = calibration_parts.mkreblock(pipeline, fccok)
		
	#
	# H(T)-OK BIT BRANCH
	#
	
	# First combine higher order bits to determine h(t)-OK
	higherbits = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, filtersok, htproduced, sciencequalitytee), "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, endianness=1234, rate=%d" % dqsr)
	higherbitstee = pipeparts.mktee(pipeline, higherbits)

	# Now calculate h(t)-OK bit
	htok = pipeparts.mkbitvectorgen(pipeline, calibration_parts.mkqueue(pipeline, higherbitstee), bit_vector = 1, threshold = 28)
	htok = pipeparts.mkcapsfilter(pipeline, htok, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, rate=%d, endianness=1234" % dqsr)
	htok = calibration_parts.mkaudiorate(pipeline, htok)
	htok = calibration_parts.mkreblock(pipeline, htok)

	#
	# HW INJECTION BITS
	#	

	hwinjcbc = pipeparts.mkgeneric(pipeline, calibration_parts.mkqueue(pipeline, odcstatevectortee), "lal_logical_undersampler", required_on = int(options.hw_inj_cbc_bitmask), status_out = 64)
	hwinjcbc = pipeparts.mkcapsfilter(pipeline, hwinjcbc, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, endianness=1234, rate=%d" % dqsr)
	hwinjcbc = calibration_parts.mkaudiorate(pipeline, hwinjcbc)
	hwinjcbc = calibration_parts.mkreblock(pipeline, hwinjcbc)

	hwinjburst = pipeparts.mkgeneric(pipeline, calibration_parts.mkqueue(pipeline, odcstatevectortee), "lal_logical_undersampler", required_on = int(options.hw_inj_burst_bitmask), status_out = 128)
	hwinjburst = pipeparts.mkcapsfilter(pipeline, hwinjburst, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, endianness=1234, rate=%d" % dqsr)
	hwinjburst = calibration_parts.mkaudiorate(pipeline, hwinjburst)
	hwinjburst = calibration_parts.mkreblock(pipeline, hwinjburst)

	hwinjdetchar = pipeparts.mkgeneric(pipeline, calibration_parts.mkqueue(pipeline, odcstatevectortee), "lal_logical_undersampler", required_on = int(options.hw_inj_detchar_bitmask), status_out = 256)
	hwinjdetchar = pipeparts.mkcapsfilter(pipeline, hwinjdetchar, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, endianness=1234, rate=%d" % dqsr)
	hwinjdetchar = calibration_parts.mkaudiorate(pipeline, hwinjdetchar)
	hwinjdetchar = calibration_parts.mkreblock(pipeline, hwinjdetchar)

	hwinjstoch = pipeparts.mkgeneric(pipeline, calibration_parts.mkqueue(pipeline, odcstatevectortee), "lal_logical_undersampler", required_on = int(options.hw_inj_stoch_bitmask), status_out = 32)
	hwinjstoch = pipeparts.mkcapsfilter(pipeline, hwinjstoch, "audio/x-raw-int, width=32, depth=32, signed=false, channels=1, endianness=1234, rate=%d" % dqsr)
	hwinjstoch = calibration_parts.mkaudiorate(pipeline, hwinjstoch)
	hwinjstoch = calibration_parts.mkreblock(pipeline, hwinjstoch)


	#
	# COMBINE ALL BITS TO MAKE GDS-CALIB_STATE_VECTOR
	#

	dqcaps = "audio/x-raw-int, rate=%d, width=32, depth=32, signed=false, channels=1, endianness=1234" % dqsr
	calibstatevector = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, nogap, higherbitstee, scienceintenttee, htok, hwinjcbc, hwinjburst, hwinjdetchar, hwinjstoch), dqcaps)
	if not options.no_kappaa:
		calibstatevector = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, calibstatevector, kaok), dqcaps)
	if not options.no_kappatst:
		calibstatevector = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, calibstatevector, ktstok), dqcaps)
	if not options.no_kappapu:
		calibstatevector = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, calibstatevector, kpuok), dqcaps)
	if not options.no_kappac:
		calibstatevector = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, calibstatevector, kcok), dqcaps)
	if not options.no_fcc:
		calibstatevector = calibration_parts.mkadder(pipeline, calibration_parts.list_srcs(pipeline, calibstatevector, fccok), dqcaps)

	calibstatevector = pipeparts.mkprogressreport(pipeline, calibstatevector, "progress_calibstatevec_%s" % instrument)
	dqtagstr = "channel-name=%s:GDS-CALIB_STATE_VECTOR, instrument=%s" % (instrument, instrument)
	calibstatevector = pipeparts.mktaginject(pipeline, calibstatevector, dqtagstr)
	calibstatevector = calibration_parts.mkaudiorate(pipeline, calibstatevector)
	calibstatevector = calibration_parts.mkreblock(pipeline, calibstatevector)

# Resample the \kappa_a channels at the specified recording sample rate and change them to single precision channels
if not options.no_kappaa:
	kaRout = pipeparts.mkaudioconvert(pipeline, kaRout)
	kaRout = pipeparts.mkcapsfilter(pipeline, kaRout, "audio/x-raw-float, channels=1, endianness=1234, width=32, rate=%d" % options.compute_factors_sr)
	kaRout = calibration_parts.resample(pipeline, kaRout, "audio/x-raw-float, channels=1, endianness=1234, width=32, rate=%d" % options.record_factors_sr)
 
	kaIout = pipeparts.mkaudioconvert(pipeline, kaIout)
	kaIout = pipeparts.mkcapsfilter(pipeline, kaIout, "audio/x-raw-float, channels=1, endianness=1234, width=32, rate=%d" % options.compute_factors_sr)
	kaIout = calibration_parts.resample(pipeline, kaIout, "audio/x-raw-float, width=32, endianness=1234, channels=1, rate=%d" % options.record_factors_sr)
# Resample the \kappa_pu channels at the specified recording sample rate and change them to single precision channels
if not options.no_kappapu:
	kpuRout = pipeparts.mkaudioconvert(pipeline, kpuRout)
	kpuRout = pipeparts.mkcapsfilter(pipeline, kpuRout, "audio/x-raw-float, channels=1, endianness=1234, width=32, rate=%d" % options.compute_factors_sr)
	kpuRout = calibration_parts.resample(pipeline, kpuRout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.record_factors_sr)
 
	kpuIout = pipeparts.mkaudioconvert(pipeline, kpuIout)
	kpuIout = pipeparts.mkcapsfilter(pipeline, kpuIout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.compute_factors_sr)
	kpuIout = calibration_parts.resample(pipeline, kpuIout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.record_factors_sr)
# Resample the \kappa_tst channels at the specified recording sample rate and change them to single precision channels
if not options.no_kappatst:
	ktstRout = pipeparts.mkaudioconvert(pipeline, ktstRout)
	ktstRout = pipeparts.mkcapsfilter(pipeline, ktstRout, "audio/x-raw-float, channels=1, endianness=1234, width=32, rate=%d" % options.compute_factors_sr)
	ktstRout = calibration_parts.resample(pipeline, ktstRout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.record_factors_sr)
 
	ktstIout = pipeparts.mkaudioconvert(pipeline, ktstIout)
	ktstIout = pipeparts.mkcapsfilter(pipeline, ktstIout, "audio/x-raw-float, channels=1, endianness=1234, width=32, rate=%d" % options.compute_factors_sr)
	ktstIout = calibration_parts.resample(pipeline, ktstIout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.record_factors_sr)
# Resample the \kappa_c channel at the specified recording sample rate and change it to a single precision channel
if not options.no_kappac:
	kcout = pipeparts.mkaudioconvert(pipeline, kcout)
	kcout = pipeparts.mkcapsfilter(pipeline, kcout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.compute_factors_sr)
	kcout = calibration_parts.resample(pipeline, kcout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.record_factors_sr)
# Resample the f_cc channel at the specified recording sample rate and change it to a single precision channel
if not options.no_fcc:
	fccout = pipeparts.mkaudioconvert(pipeline, fccout)
	fccout = pipeparts.mkcapsfilter(pipeline, fccout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.compute_factors_sr)
	fccout = calibration_parts.resample(pipeline, fccout, "audio/x-raw-float, width=32, channels=1, endianness=1234, rate=%d" % options.record_factors_sr)

# Gate the strain channel with all of the channels we want in frames
if not options.no_dq_vector:
	strain, calibstatevector = calibration_parts.gate_strain_for_output_frames(pipeline, strain, calibstatevector)
	strain, odcstatevectorout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, odcstatevectortee)
if not options.no_kappatst:
	strain, ktstRout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, ktstRout)
	strain, ktstIout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, ktstIout)
if not options.no_kappapu:
	strain, kpuRout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, kpuRout)
	strain, kpuIout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, kpuIout)
if not options.no_kappaa:
	strain, kaRout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, kaRout)
	strain, kaIout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, kaIout)
if not options.no_kappac:
	strain, kcout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, kcout)
if not options.no_fcc:
	strain, fccout = calibration_parts.gate_strain_for_output_frames(pipeline, strain, fccout)
strain = calibration_parts.mkaudiorate(pipeline, strain)
strain = calibration_parts.mkreblock(pipeline, strain)

# Gate everything with the strain channel so that no frames get written without a strain channel or any other channel already gated with the strain channel
straintee = pipeparts.mktee(pipeline, strain)
if not options.no_dq_vector:
	calibstatevector = calibration_parts.gate_other_with_strain(pipeline, calibstatevector, straintee)
	odcstatevectorout = calibration_parts.gate_other_with_strain(pipeline, odcstatevectorout, straintee)
if not options.no_kappatst:
	ktstRout = calibration_parts.gate_other_with_strain(pipeline, ktstRout, straintee)
	ktstIout = calibration_parts.gate_other_with_strain(pipeline, ktstIout, straintee)
if not options.no_kappapu:
	kpuRout = calibration_parts.gate_other_with_strain(pipeline, kpuRout, straintee)
	kpuIout = calibration_parts.gate_other_with_strain(pipeline, kpuIout, straintee)
if not options.no_kappaa:
	kaRout = calibration_parts.gate_other_with_strain(pipeline, kaRout, straintee)
	kaIout = calibration_parts.gate_other_with_strain(pipeline, kaIout, straintee)
if not options.no_kappac:
	kcout = calibration_parts.gate_other_with_strain(pipeline, kcout, straintee)
if not options.no_fcc:
	fccout = calibration_parts.gate_other_with_strain(pipeline, fccout, straintee)
	
# If the calibration factors are bad and we are NOT in science mode, set the strain channel to a gap. This is to avoid issues with reproducibility of the pipeline.  Note: If the calibration factors are bad for an extended period of time and we ARE in science mode, there will be problems.  I'm not fixing this on purpose.
# FIXME: The code below just locks up the pipeline. Think of another way to do this
strain = calibration_parts.mkqueue(pipeline, straintee)
if not options.no_dq_vector and (options.apply_kappatst or options.apply_kappac or options.apply_kappapu):
	calibstatevectortee = pipeparts.mktee(pipeline, calibstatevector)
	if options.apply_kappatst:
		ktstbadstate = pipeparts.mkstatevector(pipeline, calibration_parts.mkqueue(pipeline, calibstatevectortee), required_off = 2052)
		strain = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, strain), control = calibration_parts.mkqueue(pipeline, ktstbadstate), threshold = 1, invert_control = True)
	if options.apply_kappac:
		kcbadstate = pipeparts.mkstatevector(pipeline, calibration_parts.mkqueue(pipeline, calibstatevectortee), required_off = 4100)
		strain = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, strain), control = calibration_parts.mkqueue(pipeline, kcbadstate), threshold = 1, invert_control = True)
	if options.apply_kappapu:
		kpubadstate = pipeparts.mkstatevector(pipeline, calibration_parts.mkqueue(pipeline, calibstatevectortee), required_off = 1028)
		strain = pipeparts.mkgate(pipeline, calibration_parts.mkqueue(pipeline, strain), control = calibration_parts.mkqueue(pipeline, kpubadstate), threshold = 1, invert_control = True)

	strain = pipeparts.mkfirbank(pipeline, strain, fir_matrix = [[0,1]], time_domain = td)
	strain = calibration_parts.mkaudiorate(pipeline, strain)
	strain = calibration_parts.mkreblock(pipeline, strain)
	calibstatevector = calibration_parts.mkqueue(pipeline, calibstatevectortee)

#
# CREATE MUXER AND HOOK EVERYTHING UP TO IT
#

mux = pipeparts.mkframecppchannelmux(pipeline, None)

if options.frame_duration is not None:
        mux.set_property("frame-duration", options.frame_duration)
if options.frames_per_file is not None:
        mux.set_property("frames-per-file", options.frames_per_file)
mux.set_property("compression-scheme", options.compression_scheme)
mux.set_property("compression-level", options.compression_level)

# Link the output DQ vectors up to the muxer, if applicable
if not options.no_dq_vector:
	calibration_parts.mkqueue(pipeline, calibstatevector).get_pad("src").link(mux.get_pad("%s:%sCALIB_STATE_VECTOR%s" % (instrument, chan_prefix, chan_suffix)))
	calibration_parts.mkqueue(pipeline, odcstatevectorout).get_pad("src").link(mux.get_pad("%s:%s" % (instrument, options.dq_channel_name)))

# Link the strain branch to the muxer
strain.get_pad("src").link(mux.get_pad("%s:%sCALIB_STRAIN%s" % (instrument, chan_prefix, chan_suffix)))

# Link the real and imaginary parts of \kappa_a to the muxer
if not options.no_kappaa:
	calibration_parts.mkqueue(pipeline, kaRout).get_pad("src").link(mux.get_pad("%s:%sCALIB_KAPPA_A_REAL%s" % (instrument, chan_prefix, chan_suffix)))
	calibration_parts.mkqueue(pipeline, kaIout).get_pad("src").link(mux.get_pad("%s:%sCALIB_KAPPA_A_IMAGINARY%s" % (instrument, chan_prefix, chan_suffix)))

# Link the real and imaginary parts of \kappa_tst to the muxer
if not options.no_kappatst:
	calibration_parts.mkqueue(pipeline, ktstRout).get_pad("src").link(mux.get_pad("%s:%sCALIB_KAPPA_TST_REAL%s" % (instrument, chan_prefix, chan_suffix)))
	calibration_parts.mkqueue(pipeline, ktstIout).get_pad("src").link(mux.get_pad("%s:%sCALIB_KAPPA_TST_IMAGINARY%s" % (instrument, chan_prefix, chan_suffix)))

# Link the real and imaginary parts of \kappa_pu to the muxer
if not options.no_kappapu:
	calibration_parts.mkqueue(pipeline, kpuRout).get_pad("src").link(mux.get_pad("%s:%sCALIB_KAPPA_PU_REAL%s" % (instrument, chan_prefix, chan_suffix)))
	calibration_parts.mkqueue(pipeline, kpuIout).get_pad("src").link(mux.get_pad("%s:%sCALIB_KAPPA_PU_IMAGINARY%s" % (instrument, chan_prefix, chan_suffix)))

# Link the \kappa_c to the muxer
if not options.no_kappac:
	calibration_parts.mkqueue(pipeline, kcout).get_pad("src").link(mux.get_pad("%s:%sCALIB_KAPPA_C%s" % (instrument, chan_prefix, chan_suffix)))

# Link the f_cc to the muxer
if not options.no_fcc:
	calibration_parts.mkqueue(pipeline, fccout).get_pad("src").link(mux.get_pad("%s:%sCALIB_F_CC%s" % (instrument, chan_prefix, chan_suffix)))

if options.wings is not None:
	def clip_wings(pad, obj, (wings, start, end)):
		if isinstance(obj, gst.Buffer):
			startts = lal.LIGOTimeGPS(0, obj.timestamp)
			if startts >= (start + wings) and startts < (end - wings):
				return True
			elif startts < (start + wings) or startts >= (end - wings):
				return False
		elif isinstance(obj, gst.Event):
			return True
	mux.get_pad("src").add_data_probe(clip_wings, (lal.LIGOTimeGPS(options.wings, 0), lal.LIGOTimeGPS(int(options.gps_start_time), 0), lal.LIGOTimeGPS(int(options.gps_end_time), 0)))

def no_short_frames(pad, obj, (frame_duration)):
	if isinstance(obj, gst.Buffer):
		duration = lal.LIGOTimeGPS(0, obj.duration)
		if duration != frame_duration:
			return False
		else:
			return True
	elif isinstance(obj, gst.Event):
		return True
mux.get_pad("src").add_data_probe(no_short_frames, (lal.LIGOTimeGPS(options.frame_duration * options.frames_per_file, 0)))

mux = pipeparts.mkprogressreport(pipeline, mux, "progress_sink_%s" % instrument)

if options.write_to_shm_partition is not None:
	lvshmsink = gst.element_factory_make("gds_lvshmsink")
	lvshmsink.set_property("sync", False)
	lvshmsink.set_property("async", False)
	lvshmsink.set_property("shm-name", options.write_to_shm_partition)
	lvshmsink.set_property("num-buffers", 10)
	lvshmsink.set_property("blocksize", options.frame_size * options.frame_duration * options.frames_per_file)
	lvshmsink.set_property("buffer-mode", options.buffer_mode)
	pipeline.add(lvshmsink)
	mux.link(lvshmsink)
else:
	pipeparts.mkframecppfilesink(pipeline, mux, frame_type = options.frame_type, path = options.output_path, instrument = instrument) 

# Run pipeline

if options.write_pipeline is not None:
	pipeparts.write_dump_dot(pipeline, "%s.%s" %(options.write_pipeline, "NULL"), verbose = options.verbose)

# Seek the pipeline when necessary.  Note: The seekevent for frames is set above when command line is being parsed/sanity checked.
if options.data_source == "frames":
	datasource.do_seek(pipeline, seekevent)	
	print >>sys.stderr, "seeking GPS start and stop times ..."

if options.verbose:
	print >>sys.stderr, "setting pipeline state to playing ..."
if pipeline.set_state(gst.STATE_PLAYING) == gst.STATE_CHANGE_FAILURE:
	raise RuntimeError("pipeline failed to enter PLAYING state")
else:
	print "set to playing successfully"
if options.write_pipeline is not None:
	pipeparts.write_dump_dot(pipeline, "%s.%s" %(options.write_pipeline, "PLAYING"), verbose = options.verbose)
	
if options.verbose:
	print >>sys.stderr, "running pipeline ..."

mainloop.run()
