Request for code & logic review

[Eric S. Raymond]

The Python port of ntpdig is almost ready to land.  But there is one last
little bit of it I'm not sure I understand correctly.  I'm requesting review
of my code and assumptions.

Presently the adjustment and synch distance are calculated this way:

    def delta(self):
        return self.root_delay
    def epsilon(self):
        return self.root_dispersion
    def synchd(self):
        "Synchronization distance, estimates worst-case error in seconds"
        # This is "lambda" in NTP-speak, but that's a Python keyword 
        return abs(self.delta() - self.epsilon())
    def adjust(self):
        "Adjustment implied by this packet."
        return self.received - self.transmit_timestamp

where self.received is the receipt time of the response packet.  I do
it this way for no better reason than that some examples of Python SNTP
clients I found on the web use this formula.  But the C ntpdig code
suggests this is wrong.

Here's the C for the offset/synch-distance calculation:

--------------------------------------------------------------------------------
	/* Convert timestamps from network to host byte order */
	p_rdly = NTOHS_FP(rpkt->rootdelay);
	p_rdsp = NTOHS_FP(rpkt->rootdisp);
	NTOHL_FP(&rpkt->reftime, &p_ref);
	NTOHL_FP(&rpkt->org, &p_org);
	NTOHL_FP(&rpkt->rec, &p_rec);
	NTOHL_FP(&rpkt->xmt, &p_xmt);

	*precision = LOGTOD(rpkt->precision);

	TRACE(3, ("offset_calculation: LOGTOD(rpkt->precision): %f\n", *precision));

	/* Compute offset etc. */
	tmp = p_rec;
	L_SUB(&tmp, &p_org);
	LFPTOD(&tmp, t21);
	TVTOTS(tv_dst, &dst);
	dst.l_ui += JAN_1970;
	tmp = p_xmt;
	L_SUB(&tmp, &dst);
	LFPTOD(&tmp, t34);
	*offset = (t21 + t34) / 2.;
	delta = t21 - t34;

	// synch_distance is:
	// (peer->delay + peer->rootdelay) / 2 + peer->disp
	// + peer->rootdisp + clock_phi * (current_time - peer->update)
	// + peer->jitter;
	//
	// and peer->delay = fabs(peer->offset - p_offset) * 2;
	// and peer->offset needs history, so we're left with
	// p_offset = (t21 + t34) / 2.;
	// peer->disp = 0; (we have no history to augment this)
	// clock_phi = 15e-6; 
	// peer->jitter = LOGTOD(sys_precision); (we have no history to augment this)
	// and ntp_proto.c:set_sys_tick_precision() should get us sys_precision.
	//
	// so our answer seems to be:
	//
	// (fabs(t21 + t34) + peer->rootdelay) / 3.
	// + 0 (peer->disp)
	// + peer->rootdisp
	// + 15e-6 (clock_phi)
	// + LOGTOD(sys_precision)

	INSIST( FPTOD(p_rdly) >= 0. );
#if 1
	*synch_distance = (fabs(t21 + t34) + FPTOD(p_rdly)) / 3.
		+ 0.
		+ FPTOD(p_rdsp)
		+ 15e-6
		+ 0.	/* LOGTOD(sys_precision) when we can get it */
		;
	INSIST( *synch_distance >= 0. );
#else
	*synch_distance = (FPTOD(p_rdly) + FPTOD(p_rdsp))/2.0;
#endif
--------------------------------------------------------------------------------

That is rather nasty and in my opinion sufficient reason to shoot the C
version through the head.  Now I'm going to simplify it into pseudocode,
ignoring issues about endianness and timestamp epochs.

--------------------------------------------------------------------------------
	precision = log2(rpkt->precision);

	/* Compute offset etc. */
	t21 = rpkt->rec - rpkt->org
	t34 = now - rpkt->xmt
	offset = (t21 + t34) / 2.;

#if 1
	synch_distance = (fabs(t21 + t34) + rpkt->rdly) / 3.
		+ 0.
		+ rpkt->rdsp
		+ 15e-6
		+ 0.	/* LOGTOD(sys_precision) when we can get it */
		;
#else
	synch_distance = (pkt->rdly + pkt->p_rdsp)/2.0;
#endif
--------------------------------------------------------------------------------

I've omitted the delta calculation because it's only used in a diagnostic;
offset is what is applied to the local clock.

Note that the synch distance formula that is conditioned out is the one
I'm using.  Somebody must have thought it was correct at one time.  But if
I'm to believe the C above, the Python should read thus:

    def synchd(self):
        "Synchronization distance, estimates worst-case error in seconds"
        (self.receive_timestamp - self.origin_timestamp) \
		+ (now - self.transmit_timestamp) + self.rootdelay) / 3 \
		+ self.root_dispersion + 15e-6
    def adjust(self):
        "Adjustment implied by this packet."
        return ((self.receive_timestamp - self.origin_timestamp)
		+ (now - self.transmit_timestamp)) / 2

My requests to Daniel and anyone else who might understand these formulas:

(1) Check my logic and C translation.  A bug here would be bad.

(2) What is the authority for these formulas?  What RFC chapter and
    verse can we cite?
-- 
		<a href="http://www.catb.org/~esr/">Eric S. Raymond</a>

"Among the many misdeeds of British rule in India, history will look
upon the Act depriving a whole nation of arms as the blackest."
        -- Mohandas Gandhi, An Autobiography, pg 446