1 /*- 2 * Copyright (c) 2022 Colin Percival 3 * All rights reserved. 4 * 5 * Redistribution and use in source and binary forms, with or without 6 * modification, are permitted provided that the following conditions 7 * are met: 8 * 1. Redistributions of source code must retain the above copyright 9 * notice, this list of conditions and the following disclaimer. 10 * 2. Redistributions in binary form must reproduce the above copyright 11 * notice, this list of conditions and the following disclaimer in the 12 * documentation and/or other materials provided with the distribution. 13 * 14 * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND 15 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 16 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 17 * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE 18 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 19 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 20 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 21 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 22 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 23 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 24 * SUCH DAMAGE. 25 */ 26 27 #include <sys/cdefs.h> 28 #include <sys/param.h> 29 #include <sys/systm.h> 30 #include <sys/timetc.h> 31 #include <sys/tslog.h> 32 #include <machine/cpu.h> 33 34 /** 35 * clockcalib(clk, clkname): 36 * Return the frequency of the provided timer, as calibrated against the 37 * current best-available timecounter. 38 */ 39 uint64_t 40 clockcalib(uint64_t (*clk)(void), const char *clkname) 41 { 42 struct timecounter *tc = atomic_load_ptr(&timecounter); 43 uint64_t clk0, clk1, clk_delay, n, passes = 0; 44 uint64_t t0, t1, tadj, tlast; 45 double mu_clk = 0; 46 double mu_t = 0; 47 double va_clk = 0; 48 double va_t = 0; 49 double cva = 0; 50 double d1, d2; 51 double inv_n; 52 uint64_t freq; 53 54 TSENTER(); 55 /*- 56 * The idea here is to compute a best-fit linear regression between 57 * the clock we're calibrating and the reference clock; the slope of 58 * that line multiplied by the frequency of the reference clock gives 59 * us the frequency we're looking for. 60 * 61 * To do this, we calculate the 62 * (a) mean of the target clock measurements, 63 * (b) variance of the target clock measurements, 64 * (c) mean of the reference clock measurements, 65 * (d) variance of the reference clock measurements, and 66 * (e) covariance of the target clock and reference clock measurements 67 * on an ongoing basis, updating all five values after each new data 68 * point arrives, stopping when we're confident that we've accurately 69 * measured the target clock frequency. 70 * 71 * Given those five values, the important formulas to remember from 72 * introductory statistics are: 73 * 1. slope of regression line = covariance(x, y) / variance(x) 74 * 2. (relative uncertainty in slope)^2 = 75 * (variance(x) * variance(y) - covariance(x, y)^2) 76 * ------------------------------------------------ 77 * covariance(x, y)^2 * (N - 2) 78 * 79 * We adjust the second formula slightly, adding a term to each of 80 * the variance values to reflect the measurement quantization. 81 * 82 * Finally, we need to determine when to stop gathering data. We 83 * can't simply stop as soon as the computed uncertainty estimate 84 * is below our threshold; this would make us overconfident since it 85 * would introduce a multiple-comparisons problem (cf. sequential 86 * analysis in clinical trials). Instead, we stop with N data points 87 * if the estimated uncertainty of the first k data points meets our 88 * target for all N/2 < k <= N; this is not theoretically optimal, 89 * but in practice works well enough. 90 */ 91 92 /* 93 * Initial values for clocks; we'll subtract these off from values 94 * we measure later in order to reduce floating-point rounding errors. 95 * We keep track of an adjustment for values read from the reference 96 * timecounter, since it can wrap. 97 */ 98 clk0 = clk(); 99 t0 = tc->tc_get_timecount(tc) & tc->tc_counter_mask; 100 tadj = 0; 101 tlast = t0; 102 103 /* Loop until we give up or decide that we're calibrated. */ 104 for (n = 1; ; n++) { 105 /* Get a new data point. */ 106 clk1 = clk() - clk0; 107 t1 = tc->tc_get_timecount(tc) & tc->tc_counter_mask; 108 while (t1 + tadj < tlast) 109 tadj += (uint64_t)tc->tc_counter_mask + 1; 110 tlast = t1 + tadj; 111 t1 += tadj - t0; 112 113 /* If we spent too long, bail. */ 114 if (t1 > tc->tc_frequency) { 115 printf("Statistical %s calibration failed! " 116 "Clocks might be ticking at variable rates.\n", 117 clkname); 118 printf("Falling back to slow %s calibration.\n", 119 clkname); 120 freq = (double)(tc->tc_frequency) * clk1 / t1; 121 break; 122 } 123 124 /* Precompute to save on divisions later. */ 125 inv_n = 1.0 / n; 126 127 /* Update mean and variance of recorded TSC values. */ 128 d1 = clk1 - mu_clk; 129 mu_clk += d1 * inv_n; 130 d2 = d1 * (clk1 - mu_clk); 131 va_clk += (d2 - va_clk) * inv_n; 132 133 /* Update mean and variance of recorded time values. */ 134 d1 = t1 - mu_t; 135 mu_t += d1 * inv_n; 136 d2 = d1 * (t1 - mu_t); 137 va_t += (d2 - va_t) * inv_n; 138 139 /* Update covariance. */ 140 d2 = d1 * (clk1 - mu_clk); 141 cva += (d2 - cva) * inv_n; 142 143 /* 144 * Count low-uncertainty iterations. This is a rearrangement 145 * of "relative uncertainty < 1 PPM" avoiding division. 146 */ 147 #define TSC_PPM_UNCERTAINTY 1 148 #define TSC_UNCERTAINTY TSC_PPM_UNCERTAINTY * 0.000001 149 #define TSC_UNCERTAINTY_SQR TSC_UNCERTAINTY * TSC_UNCERTAINTY 150 if (TSC_UNCERTAINTY_SQR * (n - 2) * cva * cva > 151 (va_t + 4) * (va_clk + 4) - cva * cva) 152 passes++; 153 else 154 passes = 0; 155 156 /* Break if we're consistently certain. */ 157 if (passes * 2 > n) { 158 freq = (double)(tc->tc_frequency) * cva / va_t; 159 if (bootverbose) 160 printf("Statistical %s calibration took" 161 " %lu us and %lu data points\n", 162 clkname, (unsigned long)(t1 * 163 1000000.0 / tc->tc_frequency), 164 (unsigned long)n); 165 break; 166 } 167 168 /* 169 * Add variable delay to avoid theoretical risk of aliasing 170 * resulting from this loop synchronizing with the frequency 171 * of the reference clock. On the nth iteration, we spend 172 * O(1 / n) time here -- long enough to avoid aliasing, but 173 * short enough to be insignificant as n grows. 174 */ 175 clk_delay = clk() + (clk() - clk0) / (n * n); 176 while (clk() < clk_delay) 177 cpu_spinwait(); /* Do nothing. */ 178 } 179 TSEXIT(); 180 return (freq); 181 } 182