# Pull in the packages we need:
require(stats)
#### resmatch ####################################################
# Parameters:
#
# mfgTol - Resistor manufacturing tolerance, as a fraction. The
# default, 1/100, means we start with 1% resistors.
#
# handTol - Desired percent tolerance for the hand-matched
# resistors. Defaults to 0.1/100, or 0.1%.
#
# processFudgeFactor - Lets you simulate a more or less precise
# manufacturing process than the default, which simulates an
# ideal "Six Sigma" process, giving a few out-of-spec resistors
# per million. Higher values simulate a less precise process.
#
# resistors - Number of random resistor values to generate for the
# simulation. Defaults to a million. The lower this is, the less
# accurate the simulation, but the faster it runs. Make it too
# high and it just takes longer to run without improving
# perceptible quality.
#
# attempts - Number of matched pairs to draw from the random
# resistor set. We don't want to try matching all resistors in
# the manufacturing batch; we're simulating pulling a small
# number of matches from a larger set. Like the resistor count,
# this number is also a tradeoff between simulation accuracy
# and run time. BEWARE: If attempts is too high, the simulation
# can run out of resistors and give bogus data, because it only
# looks at each resistor once, rejecting or using each for a
# pair. If you intend to simulate a really high hand-matching
# tolerance or a really low manufacturing tolerance, you may
# need to pass a lower value for attempts or a higher value
# for resistors.
#
# nominal - Nominal resistor value, in ohms. This is only a
# scaling factor, it won't affect the sense of the results.
# The default 1Ω just keeps things simple. The only reason
# to change it is if you need the graphs to have some
# specific scale.
#
# *Color - Colors used in making PDF graphs. Defaults are for the
# http://tangentsoft.net/audio/resmatch.html graphs
resmatch <- function(
mfgTol = 1 / 100,
handTol = 0.1 / 100,
processFudgeFactor = 1,
resistors = 10 ^ 6,
attempts = resistors / 100,
nominal = 1,
barColor = "#fea300", # pumpkin
lineColor = "#501000", # milk chocolate
bgColor = "#ffffee") # cookie (getting hungry?)
{
# Common procedure for creating a PDF drawing context with the
# passed color scheme
setupPDF <- function(name) {
pdf(name)
par(bg = bgColor, fg = lineColor,
col.axis = lineColor, col.lab = lineColor,
col.main = lineColor, col.sub = lineColor)
}
# The maximum difference from the nominal value, in ohms, that
# a resistor is allowed to be and still be in spec
maxSpecError = nominal * mfgTol
# The minimum and maximum in-tolerance values. These are just
# for convenience, based on the above values.
minInTolValue = nominal - maxSpecError
maxInTolValue = nominal + maxSpecError
# Calculate the variance over which to pick random resistor
# values. The constant was hand-tweaked to give a "Six Sigma"
# result. To change the process behavior, pass a different
# processFudgeFactor; do not molest the magic constant.
var = 0.22 * maxSpecError * processFudgeFactor
# Create the resistor sample set.
cat("Running simulation...\n")
sample = sample(rnorm(resistors, nominal, var), resistors)
# Scan along data sample until we find a pair that are within
# the desired tolerance. Then, start scanning again from the
# next element past the set containing the previous match.
# This method ensures that we never reuse a resistor value, so
# we get the most out of the entropy we've created. We only do
# this for a fraction of the sample set...just enough to get a
# representative data set. You can calculate longer, but it
# probably won't change the results enough to worry about.
setStart = 1
setEnd = 2
tries = {}
while (setEnd < attempts) {
i = setStart
while (i < setEnd && setStart < setEnd) {
if (handTol >=
abs(sample[i] - sample[setEnd] / sample[i])) {
tries[[length(tries) + 1]] = setEnd - setStart + 1
setEnd = setEnd + 1
setStart = setEnd
}
i = i + 1
}
setEnd = setEnd + 1
}
# Summarize the match counts into unit bins starting from 2,
# the smallest possible match count. Then, rescale the match
# counts so that the largest value is its percentage of the
# total, and graph the results.
matches = {}
count = 2
while (count < max(tries)) {
matches[count - 1] = 0
for (m in tries) {
if (count == m) {
matches[count - 1] = matches[count-1] + 1
}
}
count = count + 1
}
# Create labels for the graphs
matches = matches / sum(matches) * 100
xlabel = 2:(length(matches) + 1)
# Open pdf device and plot matches
cat("Creating graphs: matches...\n")
setupPDF("matches.pdf")
barplot(matches, names.arg = xlabel, col = barColor,
xlab = "Match tries",
ylab = "Percent of total")
dev.off()
# Calculate cumulative probabilities of getting a match in N
# or less tries, and graph the result. The lines show where
# the probabilities of a good match are greater than 90%, 95%,
# and 99%.
probaccum = cumsum(matches)
# Open PDF device and plot the cumulative distribution and
# important percentage lines
cat(" probability accumulation...\n")
setupPDF("probaccum.pdf")
barplot(probaccum, names.arg = xlabel, col = barColor,
xlab = "Match tries",
ylab = "% chance of good match")
abline(a = 90, b = 0)
abline(a = 95, b = 0)
abline(a = 99, b = 0)
dev.off()
# Graph the data frequency as a visual check that the above
# parameters are generating a reasonable data set. Also, it's
# pretty.
cat(" data quality check histogram...\n")
setupPDF("histogram.pdf")
hist(sample, freq = FALSE, col = barColor)
lines(density(sample))
dev.off()
# Check how many outliers we have as a further check on the
# data's quality. The first result is the number of below-spec
# resistors, the second is the number of those above spec. Six
# sigma means 3.4 outliers per million, but as long as both
# values are in the single digits, it's close enough. I don't
# like these values to get too low: if it says 0, how can you
# tell how much better than six sigma you are?
belowSpec = 0
aboveSpec = 0
for (i in sample) {
if ((i - minInTolValue) <= 0) {
belowSpec = belowSpec+1
}
if ((i - maxInTolValue) >= 0) {
aboveSpec = aboveSpec+1
}
}
cat("Number of resistors below spec: ", belowSpec, "\n")
cat("Number of resistors above spec: ", aboveSpec, "\n")
}
resmatch()