diff --git a/gstlal-inspiral/python/stats/inspiral_extrinsics.py b/gstlal-inspiral/python/stats/inspiral_extrinsics.py
index 23804ead6fb7c5c70848a62d84bdc227542a75bc..7fddde91860481fe6d08c815f6f2644d3b8414d4 100644
--- a/gstlal-inspiral/python/stats/inspiral_extrinsics.py
+++ b/gstlal-inspiral/python/stats/inspiral_extrinsics.py
@@ -1144,9 +1144,15 @@ class TimePhaseSNR(object):
 	locations = {"H1":lal.CachedDetectors[lal.LHO_4K_DETECTOR].location, "L1":lal.CachedDetectors[lal.LLO_4K_DETECTOR].location, "V1":lal.CachedDetectors[lal.VIRGO_DETECTOR].location}#, "K1":lal.CachedDetectors[lal.KAGRA_DETECTOR].location}
 	numchunks = 20
 
-	# FIXME compute this more reliably or expose it as a property
-	# or something
-	sigma = {"time": 0.001, "phase": numpy.pi / 6, "deff": 0.2}
+	# These are generated by running
+	# ./gstlal_inspiral_compute_dtdphideff_cov_matrix --psd-xml 2019-05-09-H1L1V1psd_new.xml.gz --H-snr 5.0 --L-snr 7.0 --V-snr 2.25
+	# NOTE NOTE NOTE: You cannot just change these without fully rebuilding
+	# the trees and the marginalized distributions.
+	transtt = {frozenset(("H1", "L1")): 5.19313416e+04, frozenset(("H1", "V1")): 1.39663777e+04, frozenset(("L1", "V1")): 3.98320074e+04}
+	transpp = {frozenset(("H1", "L1")): 3.50152221e+00, frozenset(("H1", "V1")): 1.56080114e+00, frozenset(("L1", "V1")): 1.63248925e+00}
+	transtp = {frozenset(("H1", "L1")): -6.51410977e+01, frozenset(("H1", "V1")): -2.07581116e+01, frozenset(("L1", "V1")): -6.05964763e+01}
+	transpt = {frozenset(("H1", "L1")): 0.0, frozenset(("H1", "V1")): 0.0, frozenset(("L1", "V1")): 0.0}
+	transdd = {frozenset(("H1", "L1")): 5.0 , frozenset(("H1", "V1")): 5.0, frozenset(("L1", "V1")): 5.0}
 
 	def __init__(self, tree_data = None, margsky = None, verbose = False, margstart = 0, margstop = None):
 		"""
@@ -1358,7 +1364,7 @@ class TimePhaseSNR(object):
 		"""
 		Given a dictionary of time, phase and deff, which could be
 		lists of values, pack the delta t delta phi and eff distance
-		ratios divided by the values in self.sigma into an output array according to
+		ratios transformed by the covariance matrix according to
 		the rules provided by slices.
 
 		>>> TimePhaseSNR.dtdphideffpoints({"H1":0, "L1":-.001, "V1":.001}, {"H1":0, "L1":0, "V1":1}, {"H1":1, "L1":3, "V1":4}, TimePhaseSNR.slices)
@@ -1379,10 +1385,16 @@ class TimePhaseSNR(object):
 
 		for ifos, slc in slices.items():
 			ifo1, ifo2 = ifos
-			out[:,slc[0]] = (time[ifo1] - time[ifo2]) / self.sigma["time"]
-			out[:,slc[1]] = (phase[ifo1] - phase[ifo2]) / self.sigma["phase"]
+			dt = (time[ifo1] - time[ifo2])
+			dphi = (phase[ifo1] - phase[ifo2])
+			# FIXME precompute?
+			dtdphivec = numpy.array([dt, dphi])
+			coordtransmat = numpy.array([[self.transtt[frozenset((ifo1, ifo2))], self.transtp[frozenset((ifo1, ifo2))]],[self.transpt[frozenset((ifo1, ifo2))], self.transpp[frozenset((ifo1, ifo2))]]])
+			out[:,slc[0]], out[:,slc[1]] = numpy.dot(coordtransmat, dtdphivec)
+			#out[:,slc[0]] = self.transtt[frozenset((ifo1, ifo2))] * dt + dphi * self.transpt[frozenset((ifo1, ifo2))]
+			#out[:,slc[1]] = self.transpp[frozenset((ifo1, ifo2))] * dt + self.transtp[frozenset((ifo1, ifo2))] * dphi
 			# FIXME should this be the ratio - 1 or without the -1 ???
-			out[:,slc[2]] = (deff[ifo1] / deff[ifo2] - 1) / self.sigma["deff"]
+			out[:,slc[2]] = (deff[ifo1] / deff[ifo2] - 1) * self.transdd[frozenset((ifo1, ifo2))]
 
 		return out