Thursday, April 9, 2015

Half Normal Probability Paper


The code below can be used to generate half normal probability paper.  The erfinv function is a bit inaccurate (2 sig figs) so it could be improved with a better approximation.  Update: I modified this with the code from Wikipedia to give better accuracy.   Create a blank file first.  Mine was 1500 by 2000.  The code is Python.  I had to install Pillow for this.

from PIL import Image,ImageDraw,ImageFont
from math import *

def halfnormal(z,s):
  return erf(z/sqrt(2)/s)

def erfinv(z):
  a = .147
  x=sqrt(sqrt((2/pi/a+log(1-z**2)/2)**2-log(1-z**2)/a)-(2/pi/a+log(1-z**2)/2))
  #print z,z-erf(x)
  return x
  
def invhalfnorm(z,s):
  return erfinv(z)*sqrt(2)
  
def halfnormalgraph():
  im = Image.open("blank_graph.png")
  draw = ImageDraw.Draw(im)
  font1 = ImageFont.truetype("arial.ttf", 16)
  font2 = ImageFont.truetype("arial.ttf", 28)
  j=30.0
  k=700.0
  xmax=46
  xshift=30
  gridfill = (180,180,250)
  
  draw.text((550,10), "Half Normal Probability Paper", font=font2, fill=256)
  for i in range(xmax):
    x1=j*i+j*2+xshift
    y1=im.size[1]-j*2
    y2=j*2
    draw.line((x1, y1, x1, y2), fill=gridfill)
  x1=j*2+xshift
  x2=j*xmax+j+xshift
  for i in range(10):
    iv=i/20.0+0
    y1=im.size[1]-invhalfnorm(iv,1.0)*k-j*2
    draw.line((x1, y1, x2, y1 ), fill=gridfill)
    draw.text((2+xshift,y1-8), "{:10.0f}".format(iv*100), font=font1, fill=256)
  for i in range(20):
    iv=i/50.0+.5
    y1=im.size[1]-invhalfnorm(iv,1.0)*k-j*2
    draw.line((x1, y1, x2, y1 ), fill=gridfill)
    draw.text((2+xshift,y1-8), "{:10.0f}".format(iv*100), font=font1, fill=256)
  for i in range(10):
    iv=i/100.0+.90
    y1=im.size[1]-invhalfnorm(iv,1.0)*k-j*2
    draw.line((x1, y1, x2, y1 ), fill=gridfill)
    draw.text((2+xshift,y1-8), "{:10.0f}".format(iv*100), font=font1, fill=256)
  im.save("halfnormgraph.png")

Here is another program that makes the desired number of divisions and an example of the output.  This would be useful for plotting full factorial designed experiments, DOE, with 2 levels and 3 factors.  Included are sheets for 7, 15, and 31 divisions which are the ones most useful for full factorial experiments.





def hng2(bins=7):
  im = Image.open("blank_graph.png")
  draw = ImageDraw.Draw(im)
  font1 = ImageFont.truetype("arial.ttf", 16)
  font2 = ImageFont.truetype("arial.ttf", 28)
  j=30.0
  k=700.0
  xmax=46
  xshift=30
  gridfill = (180,180,250)
  
  draw.text((550,10), "Half Normal Probability Paper", font=font2, fill=256)
  for i in range(xmax):
    x1=j*i+j*2+xshift
    y1=im.size[1]-j*2
    y2=j*2
    draw.line((x1, y1, x1, y2), fill=gridfill)
  x1=j*2+xshift
  x2=j*xmax+j+xshift
  iv=0.0
  y1=im.size[1]-invhalfnorm(iv,1.0)*k-j*2
  draw.line((x1, y1, x2, y1 ), fill=gridfill)
  draw.text((-2+xshift,y1-8), "{:10.2f}".format(iv*100), font=font1, fill=256)
  for i in range(bins):
    iv=1.0*i/bins+1.0/2.0/bins
    y1=im.size[1]-invhalfnorm(iv,1.0)*k-j*2
    draw.line((x1, y1, x2, y1 ), fill=gridfill)
    draw.text((-2+xshift,y1-8), "{:10.2f}".format(iv*100), font=font1, fill=256)
  im.save("hng2.png")

hng2(7)

No comments:

Post a Comment