diff --git a/.gitignore b/.gitignore index b25c15b..518a1c6 100644 --- a/.gitignore +++ b/.gitignore @@ -1 +1,2 @@ *~ +*.log diff --git a/infocalc.py b/infocalc.py new file mode 100755 index 0000000..6770d3b --- /dev/null +++ b/infocalc.py @@ -0,0 +1,113 @@ +#!/usr/bin/env python +## Axiomatic Design Information Calculator (and plotter) +## Author: Joseph Timothy foley +## Start Date: 2026-02-27 +## Input: data in csv file +## Output: information calculation and PDF for report/presentation +import os +import logging +import argparse +from pathlib import PurePath##https://docs.python.org/3/library/pathlib.html#module-pathlib +import numpy as np +import matplotlib.pyplot as plt +from scipy.stats import norm +import pandas as pd + +#Main program loop +print("""Axiomatic Design Information Calculator by Joseph. T. Foley +From https://gitea.cs.ru.is/AxiomaticDesign/adcalc/""") +parser = argparse.ArgumentParser( + description="Axiomatic Design Information Calculator.") +parser.add_argument('csvfile', + help="CSV file with data and headers") +parser.add_argument('column', + help='Which column header to take data from') +parser.add_argument('minvalue', type=float, + help='Tolerance low limit') +parser.add_argument('maxvalue', type=float, + help='Tolerance high limit') +parser.add_argument('--normalizey', action="store_true", + help='Set y-axis to normalized probability density') +parser.add_argument('--log', default="INFO", + help='Console log level: Number or DEBUG, INFO, WARNING, ERROR') +parser.add_argument('--graphinfo', + help='Put information on the PDF graph') +args = parser.parse_args() + +## Set up logging +numeric_level = getattr(logging, args.log.upper(), None) +if not isinstance(numeric_level, int): + raise ValueError(f'Invalid log level: {args.log}') +#print(f"Log level: {numeric_level}") +logger = logging.getLogger("app") +logger.setLevel(numeric_level) +# log everything to file +logpath = os.path.splitext(args.csvfile)[0]+".log" +fh = logging.FileHandler(logpath) +fh.setLevel(logging.DEBUG) +# log to console +ch = logging.StreamHandler() +ch.setLevel(numeric_level) +# create formatter and add to handlers +consoleformatter = logging.Formatter('%(message)s') +ch.setFormatter(consoleformatter) +spamformatter = logging.Formatter('%(asctime)s %(name)s[%(levelname)s] %(message)s') +fh.setFormatter(spamformatter) +# add the handlers to logger +logger.addHandler(ch) +logger.addHandler(fh) + +logger.info("Creating infocalc log file %s", logpath) + +# filename pre-processing for output +inpath = PurePath(args.csvfile) +print(f"Input: {inpath}") + +# grab the data and process +data = np.array(pd.read_csv(inpath)[args.column]) +lowerbound = args.minvalue +upperbound = args.maxvalue +logger.debug(f"data:{data}, lower:{lowerbound}, upper:{upperbound}") + +mean = data.mean() +stddev = data.std(ddof=1) +# Delta Degrees of Freedom: ddof=0 for population, ddof=1 for sample std dev +prob = norm.cdf(upperbound, mean, stddev) - norm.cdf(lowerbound, mean, stddev) +#print("probability: %f", prob) +info = -np.emath.log2(prob) +#print("information content: %f bits", info) +## place text on plot: https://matplotlib.org/3.3.4/gallery/recipes/placing_text_boxes.html +fig, ax = plt.subplots() +textstr = '\n'.join(( + r'$n=%d$' % (len(data)), + r'$\mu=%.2f$' % (mean, ), + r'$\sigma=%.2f$' % (stddev, ), + r'$P=%.2f$' % (prob, ), + r'$I=%.2f$ bits' % (info, ))) +# these are matplotlib.patch.Patch properties +props = dict(boxstyle='round', facecolor='wheat', alpha=0.5) + +# place a text box in upper left in axes coords +ax.text(0.05, 0.95, textstr, transform=ax.transAxes, fontsize=14, + verticalalignment='top', bbox=props) + +x = np.linspace(mean-3*stddev, mean+3*stddev, 500) +y = norm.pdf(x, loc=mean, scale=stddev) +if args.normalizey: + y = y * stddev#rescale back to unity area +plt.axvline(x=mean, color="green", linestyle="dashed", label="mean") +plt.axvline(lowerbound, color="red") +plt.axvline(upperbound, color="red") + +plt.plot(x, y, 'b-', label='Normal distribution') +#yt = scipy.stats.t.pdf(x, len(data)-1, mean, stddev) +#plt.plot(x, yt, 'g-', label='T Distribution') +coloredregion = (x >= lowerbound) & ( x <= upperbound ) #select x values +plt.fill_between(x, 0, y, where=coloredregion, color="grey", alpha=0.5, label="Design range") +plt.xlabel('X') +plt.ylabel('Probability density') +plt.legend() +plt.grid(True) + +top = plt.ylim()[1] +plt.show()