%% $Id: pst-fit.tex 673 2012-04-01 09:50:48Z herbert $ %% %% This is file `pst-fit.tex', %% %% IMPORTANT NOTICE: %% %% Package `pst-fit.tex' %% %% Buddy Ledger %% Herbert Voss %% %% This program can be redistributed and/or modified under the terms %% of the LaTeX Project Public License Distributed from CTAN archives %% in directory macros/latex/base/lppl.txt. %% %% DESCRIPTION: %% `pst-fit' is a PSTricks package for curve fitting %% \csname PSTfitLoaded\endcsname \let\PSTfitLoaded\endinput % Requires PSTricks, pst-plot \ifx\PSTricksLoaded\endinput\else \input pstricks.tex\fi \ifx\PSTricksAddLoaded\endinput\else\input pstricks-add.tex\fi \ifx\PSTXKeyLoaded\endinput\else \input pst-xkey \fi % \def\fileversion{0.02} \def\filedate{2017/08/24} \message{`PST-fit' v\fileversion, \filedate\space (BL,HV)} % \edef\PstAtCode{\the\catcode`\@} \catcode`\@=11\relax % \pst@addfams{pst-fit} % \define@key[psset]{pst-fit}{symbolFont}{\def\psk@symbolFont{/#1 }}% \psset[pst-fit]{symbolFont=StandardSymL} % \def\beginplot@ErrorLine{\begin@OpenObj} \def\endplot@ErrorLine{\psErrorLine@ii} \def\psErrorLine@ii{% \addto@pscode{ \pst@cp % current point \psline@iii % arc and lineto type \tx@ErrorLine }% \end@OpenObj% } % \def\tx@ErrorLine{ErrorLine } \pst@def{ErrorLine}<{ /copypt { dup 3 2 roll dup 4 1 roll exch } def %copy point NArray n 0 eq not { n { \psk@relyerr\space 0 eq not {copypt copypt dup \psk@relyerr\space mul add ArrowA dup \psk@relyerr\space mul sub Lineto copypt dup \psk@relyerr\space mul sub ArrowB L} if \psk@relxerr\space 0 eq not {copypt copypt exch dup \psk@relxerr\space mul add exch ArrowA exch dup \psk@relxerr\space mul sub exch Lineto copypt exch dup \psk@relxerr\space mul sub exch ArrowB L pop pop } { pop pop } ifelse } repeat } if }> \define@boolkey[psset]{pst-fit}[Pst@]{ScyBase}[true]{}% \define@boolkey[psset]{pst-fit}[Pst@]{ScxBase}[true]{}% \psset{ScyBase=true,ScxBase=true} \def\pstScalePoints(#1,#2)#3#4{% % xScale | yScale | xOperator | yOperator % the operators can be any Postscript code \pst@def{ScalePoints}<% /y ED /x ED /yOper { #4\space y mul #2\space mul } def /xOper { #3\space x mul #1\space mul } def counttomark dup dup cvi eq not { exch pop } if /m exch def /n m 2 div cvi def n { #4\space y mul #2\space mul m 1 roll #3\space x mul #1\space mul m 1 roll /m m 2 sub def } repeat>% } \def\psLineToXAxis@ii{% \addto@pscode{\pst@cp \psline@iii \tx@LineToXAxis}% \end@OpenObj} \def\tx@LineToXAxis{LineToXAxis } \pst@def{LineToXAxis}<{% NArray n 0 eq not { %n 1 eq { 0 0 /n 2 def } if %ArrowA %/n n 2 sub def %2 copy pop 0 \ifPst@ScyBase yOper \fi moveto Lineto n { 2 copy 2 copy pop 0 \ifPst@ScyBase yOper \fi 4 2 roll ArrowA Lineto 2 copy pop 0 \ifPst@ScyBase yOper \fi ArrowB L pop pop } repeat %2 copy moveto pop 0 \ifPst@ScyBase yOper \fi %L %pop pop } if}> \def\psLineToYAxis@ii{% \addto@pscode{\pst@cp \psline@iii \tx@LineToYAxis}% \end@OpenObj} \def\tx@LineToYAxis{LineToYAxis } \pst@def{LineToYAxis}<{ NArray n 0 eq not { %n 1 eq { 0 0 /n 2 def } if %ArrowA %/n n 2 sub def %CP 2 copy moveto exch pop 0 \ifPst@ScxBase xOper \fi exch Lineto n { 2 copy 2 copy exch pop 0 \ifPst@ScxBase xOper \fi exch 4 2 roll ArrowA Lineto 2 copy exch pop 0 \ifPst@ScxBase xOper \fi exch ArrowB L pop pop } repeat %2 copy moveto exch pop 0 \ifPst@ScxBase xOper \fi exch %L %pop pop } if}> % %%%%%%PREPAREPOINTS \define@key[psset]{pst-fit}{plotNoTwo}{\def\psk@plotNoTwo{#1}}% \define@key[psset]{pst-fit}{plotNoTwoFunc}{\def\psk@plotNoTwoFunc{#1}} %Must operate on two stack elements and return one \psset{plotNoTwo=0,plotNoTwoFunc=dup mul exch dup mul add sqrt} %\def\tx@PreparePoints{PreparePoints } \pst@def{PreparePoints}<{% counttomark /m exch def \ifPsk@xyValues\else % we have only y values /mm m def /M m 1 add def m { mm exch M 2 roll /M M 1 add def /mm mm 1 sub def } repeat /m m dup add def \fi \ifPst@ChangeOrder /m0 m def m \psk@plotNoMax\space 1 add div 1 sub cvi { m0 \psk@plotNoMax\space 1 add roll /m0 m0 \psk@plotNoMax\space 1 add sub def } repeat \fi /n m \psk@plotNoMax\space 1 add div cvi def \psk@plotNoMax\space 1 gt {% multiple data files? n { \psk@plotNoTwo\space 0 gt { \psk@plotNoMax\space copy \psk@plotNoMax\space \psk@plotNoTwo\space 1 sub neg roll % x yNo y yNoTwo y [yNo y yNoTwo y] \psk@plotNoMax\space 1 sub { pop } repeat % x yNo y yNoTwo y yNoTwo \psk@plotNoMax\space 1 add 1 roll } if \psk@plotNoMax\space \psk@plotNo\space 1 sub neg roll % x yNo y y y ... \psk@plotNoMax\space 1 sub { pop } repeat % x [yNoTwo] yNo \psk@plotNoTwo\space 0 gt { \psk@plotNoTwoFunc\space } if %x f(yNo,yNoTwo) /m m \psk@plotNoMax\space 1 sub sub def m 2 roll } repeat } if % no multiple data files % counttomark /m exch def % /n m 2 div cvi def /xMax -99999 def /yMax -99999 def /xP 0 def /yP 0 def m copy n { /y exch def /x exch def xMax x lt { /xMax x def } if yMax y lt {/yMax y def } if xP x gt { /xP x def } if yP y gt { /yP y def } if } repeat % m 2 roll \psk@xStep\space 0 gt \psk@yStep\space 0 gt or (\psk@xStart) length 0 gt or (\psk@yStart) length 0 gt or (\psk@xEnd) length 0 gt or (\psk@yEnd) length 0 gt or { % (\psk@xStart) length 0 gt {\psk@xStart\space }{ xP } ifelse /xStart exch def (\psk@yStart) length 0 gt {\psk@yStart\space }{ yP } ifelse /yStart exch def (\psk@xEnd) length 0 gt { \psk@xEnd\space }{ xMax } ifelse /xEnd exch def (\psk@yEnd) length 0 gt { \psk@yEnd\space }{ yMax } ifelse /yEnd exch def n { m -2 roll 2 copy /yVal exch def /xVal exch def xVal xP ge yVal yP ge and xVal xEnd le and yVal yEnd le and xVal xStart ge and yVal yStart ge and { /xP xP \psk@xStep\space add def /yP yP \psk@yStep\space add def }{% pop pop /m m 2 sub def } ifelse } repeat }{% /ncount 0 def (\psk@nEnd) length 0 gt { \psk@nEnd\space }{ m } ifelse /nEnd exch def n { m -2 roll \psk@nStep\space 1 gt { ncount \psk@nStart\space sub \psk@nStep\space mod 0 eq }{ true } ifelse ncount nEnd le and ncount \psk@nStart\space ge and not { pop pop /m m 2 sub def } if /ncount ncount 1 add def } repeat } ifelse }> % \def\beginplot@PrintCoor{\begin@SpecialObj} \def\endplot@PrintCoor{% \psPrintCoor@ii %\psk@fillstyle\ifpsshadow\pst@closedshadow\fi% %\pst@stroke \end@SpecialObj% } \define@key[psset]{pst-fit}{relxerr}{\def\psk@relxerr{#1}} \define@key[psset]{pst-fit}{relyerr}{\def\psk@relyerr{#1}} \define@key[psset]{pst-fit}{yShift}{\def\psk@yShift{#1}} \define@boolkey[psset]{pst-fit}[Pst@]{science}[true]{}% \psset[pst-fit]{science=false,relxerr=0.1,relyerr=0.1,yShift=0} %Following is new pst-func.tex key which conflicts with mine. %Must specify [pst-fit] in psset to correct. %\define@boolkey[psset]{pst-add}[Pst@]{science}[true]{% % \ifPst@science\def\psk@Scin{true }\else\def\psk@Scin{false }\fi} \def\psPrintCoor@ii{% \addto@pscode{% false \tx@NArray \pst@number\psxunit \pst@number\psyunit \psPrintCoor@iii }} \def\psPrintCoor@iii{% /mfont { \psk@PSfont\space findfont \psk@fontscale\space scalefont setfont } bind def /mfontexp { \psk@PSfont\space findfont \psk@fontscale 1.2 div scalefont setfont } bind def /s1 { \psk@symbolFont\space findfont \psk@fontscale\space scalefont setfont } bind def /ttxspc \psk@fontscale\space 4 div def /txspc ttxspc 4 div def %pop /ysc ED /xsc ED n 0 eq not {% n {% ysc div /Yval ED xsc div /Xval ED \ifPst@science %X Coordinate Xval 0 gt { Xval log floor }{ Xval 0 lt { Xval -1 mul log floor }{ 0 } ifelse } ifelse cvi /xcoorexp ED Xval 10 xcoorexp exp div \psk@decimals\space -1 gt {% 10 \psk@decimals\space exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /xcoor ED %Y Coordinate Yval 0 gt { Yval log floor }{ Yval 0 lt { Yval -1 mul log floor }{ 0 } ifelse } ifelse cvi /ycoorexp ED Yval 10 ycoorexp exp div \psk@decimals\space -1 gt {% 10 \psk@decimals\space exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /ycoor exch def %Prep Strings xcoorexp \psk@valuewidth\space string cvs /xcoorexp ED ycoorexp \psk@valuewidth\space string cvs /ycoorexp ED xcoor \psk@valuewidth\space string cvs /xcoor ED ycoor \psk@valuewidth\space string cvs /ycoor ED %X Error \psk@relxerr\space Xval mul /delx ED \psk@relxerr\space 0 ne {% delx 0 gt { delx log floor }{ delx 0 lt { delx -1 mul log floor }{ 0 } ifelse } ifelse cvi /delxexp ED delx 10 delxexp exp div dup mul sqrt \psk@decimals\space -1 gt {% 10 \psk@decimals\space exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /delx ED %Prep Strings delxexp \psk@valuewidth\space string cvs /delxexp ED delx \psk@valuewidth\space string cvs /delx ED } if %Y Error \psk@relyerr\space Yval mul /dely ED \psk@relyerr\space 0 ne {% dely 0 gt { dely log floor }{ dely 0 lt { dely -1 mul log floor }{ 0 } ifelse } ifelse cvi /delyexp ED dely 10 delyexp exp div dup mul sqrt \psk@decimals\space -1 gt {% 10 \psk@decimals\space exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /dely ED %Prep Strings delyexp \psk@valuewidth\space string cvs /delyexp ED dely \psk@valuewidth\space string cvs /dely ED } if \else Xval \psk@decimals\space -1 gt { 10 \psk@decimals\space exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq {cvi} if /xcoor ED Yval \psk@decimals\space -1 gt { 10 \psk@decimals\space exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq {cvi} if /ycoor ED xcoor \psk@valuewidth\space string cvs /xcoor ED ycoor \psk@valuewidth\space string cvs /ycoor ED % \psk@relxerr\space Xval mul dup mul sqrt \psk@decimals\space -1 gt { 10 \psk@decimals\space exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq {cvi} if /delx ED \psk@relyerr\space Yval mul dup mul sqrt \psk@decimals\space -1 gt { 10 \psk@decimals\space exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq {cvi} if /dely ED delx \psk@valuewidth\space string cvs /delx ED dely \psk@valuewidth\space string cvs /dely ED \fi % newpath Xval xsc mul %\psk@xShift\space xsc mul add Yval ysc mul %\psk@yShift\space ysc mul add moveto \Pst@Debug\space 0 gt { mfont Xval \psk@valuewidth\space string cvs show (,) show Yval \psk@valuewidth\space string cvs show } {s1 (\string\050) show mfont xcoor show \ifPst@science mfont (x) show mfont (10) show 0 ttxspc rmoveto mfontexp xcoorexp show \fi \psk@relxerr\space 0 ne { \ifPst@science 0 ttxspc neg rmoveto s1 (\string\261) show mfont delx show mfont (x) show mfont (10) show 0 ttxspc rmoveto mfontexp delxexp show 0 ttxspc neg rmoveto s1 (,) show \else s1 (\string\261) show mfont delx show s1 (,) show \fi } { \ifPst@science 0 ttxspc neg rmoveto s1 (,) show \else s1 (,) show \fi } ifelse %start on ycoor mfont ycoor show \ifPst@science mfont (x) show mfont (10) show 0 ttxspc rmoveto mfontexp ycoorexp show \fi \psk@relyerr\space 0 ne { \ifPst@science 0 ttxspc neg rmoveto s1 (\string\261) show mfont dely show mfont (x) show mfont (10) show 0 ttxspc rmoveto mfontexp delyexp show 0 ttxspc neg rmoveto s1 (\string\051) show \else s1 (\string\261) show mfont dely show s1 (\string\051) show \fi } { s1 (\string\051) show } ifelse } ifelse } repeat } if } %%%POLY FIT \define@key[psset]{pst-fit}{EqPos}[{}]{\def\psk@EqPos{#1}} \define@key[psset]{pst-fit}{MaPos}[{}]{\def\psk@MaPos{#1}}%added Feb3, 2012 Matrix out \define@key[psset]{pst-fit}{MaScale}[{}]{\def\psk@MaScale{#1}}%added Feb3, 2012 Matrix out \define@boolkey[psset]{pst-fit}[Pst@]{ShowEq}[true]{}% \psset{EqPos=,ShowEq=true,MaPos=0 0,MaScale=1} \define@key[psset]{pst-fit}{PolyOrder}[{}]{\def\psk@PolyOrder{#1}} \define@boolkey[psset]{pst-fit}[Pst@]{ReduceOrder}[true]{%these percents are critical removal results in shift to in x-dir \ifPst@ReduceOrder% \def\psk@ReduceOrder{true}\else\def\psk@ReduceOrder{false}\fi}% \define@boolkey[psset]{pst-fit}[Pst@]{PowerFit}[true]{%these percents are critical removal results in shift to in x-dir \ifPst@PowerFit% \def\psk@PowerFit{true}\def\psk@LogEFit{false}% \def\psk@LogTFit{false}\def\psk@ExpFit{false}% \def\psk@GaussFit{false}\def\psk@CustomFit{false}\def\psk@RecipFit{false}% \else\def\psk@PowerFit{false}\fi}% \define@boolkey[psset]{pst-fit}[Pst@]{LogEFit}[true]{% \ifPst@LogEFit% \def\psk@PowerFit{false}\def\psk@LogEFit{true}% \def\psk@LogTFit{false}\def\psk@ExpFit{false}% \def\psk@GaussFit{false}\def\psk@CustomFit{false}\def\psk@RecipFit{false}% \else\def\psk@LogEFit{false}\fi}% \define@boolkey[psset]{pst-fit}[Pst@]{LogTFit}[true]{% \ifPst@LogTFit% \def\psk@PowerFit{false}\def\psk@LogEFit{false}% \def\psk@LogTFit{true}\def\psk@ExpFit{false}% \def\psk@GaussFit{false}\def\psk@CustomFit{false}\def\psk@RecipFit{false}% \else\def\psk@LogTFit{false}\fi}% \define@boolkey[psset]{pst-fit}[Pst@]{ExpFit}[true]{% \ifPst@ExpFit% \def\psk@PowerFit{false}\def\psk@LogEFit{false}% \def\psk@LogTFit{false}\def\psk@ExpFit{true}% \def\psk@GaussFit{false}\def\psk@CustomFit{false}\def\psk@RecipFit{false}% \else\def\psk@ExpFit{false}\fi}% \define@boolkey[psset]{pst-fit}[Pst@]{GaussFit}[true]{% \ifPst@GaussFit% \def\psk@PowerFit{false}\def\psk@LogEFit{false}% \def\psk@LogTFit{false}\def\psk@ExpFit{false}% \def\psk@GaussFit{true}\def\psk@CustomFit{false}\def\psk@RecipFit{false}% \else\def\psk@GaussFit{false}\fi}% \define@boolkey[psset]{pst-fit}[Pst@]{RecipFit}[true]{% \ifPst@RecipFit% \def\psk@PowerFit{false}\def\psk@LogEFit{false}% \def\psk@LogTFit{false}\def\psk@ExpFit{false}% \def\psk@GaussFit{false}\def\psk@CustomFit{false}\def\psk@RecipFit{true}% \else\def\psk@RecipFit{false}\fi}% \define@boolkey[psset]{pst-fit}[Pst@]{CustomFit}[true]{% \ifPst@CustomFit% \def\psk@PowerFit{false}\def\psk@LogEFit{false}% \def\psk@LogTFit{false}\def\psk@ExpFit{false}% \def\psk@GaussFit{false}\def\psk@CustomFit{true}\def\psk@RecipFit{false}% \else\def\psk@CustomFit{false}\fi}% \define@key[psset]{pst-fit}{FXtrans}[{}]{\def\psk@FXtrans{#1}} \define@key[psset]{pst-fit}{FYtrans}[{}]{\def\psk@FYtrans{#1}} \define@key[psset]{pst-fit}{RYtrans}[{}]{\def\psk@RYtrans{#1}} \define@key[psset]{pst-fit}{Yint}[{}]{\def\psk@Yint{#1}} \define@key[psset]{pst-fit}{CheckZeroX}[{}]{\def\psk@CheckZeroX{#1}} \define@key[psset]{pst-fit}{CheckZeroY}[{}]{\def\psk@CheckZeroY{#1}} \psset{PolyOrder=1,ReduceOrder=false,PowerFit=false,LogEFit=false,LogTFit=false,% ExpFit=false,GaussFit=false,RecipFit=false,CustomFit=false,% FXtrans=,FYtrans=,RYtrans=,Yint=,CheckZeroX=,CheckZeroY=} \def\beginplot@GLLSR{\begin@SpecialObj} \def\endplot@GLLSR{% \psGLLSR@ii\psk@fillstyle\ifpsshadow\pst@closedshadow\fi% \pst@stroke \end@SpecialObj% } \def\psGLLSR@ii{\addto@pscode{false \tx@NArray \psGLLSR@iii }} \def\psGLLSR@iii{% %General Fit Model: y = a0 + a1*x + a2*x^2 + ... + an*x^n %Define Fonts /mfont { \psk@PSfont\space findfont \psk@fontscale\space scalefont setfont } bind def /mfontbig { \psk@PSfont\space findfont \psk@fontscale\space 2 mul scalefont setfont } bind def /mfontexp { \psk@PSfont\space findfont \psk@fontscale\space 2 div scalefont setfont } bind def mfont %Default parameters (\psk@PolyOrder) length 0 gt { \psk@PolyOrder\space /order ED }{ /order 1 def } ifelse /PolyEq true def (\psk@Yint) length 0 gt {\psk@Yint\space /Yint ED /setYint true def /ReduceOrder true def } { /Yint 0 def /setYint false def /ReduceOrder false def} ifelse /CheckZeroX false def /CheckZeroY false def /FXtrans {1 mul} def /FYtrans { Yint sub } def /RYtrans {1 mul} def /xPr {(x) show} def %Special Fit Functions (Linear Transforms) \psk@PowerFit\space { /order 1 def /CheckZeroX true def /CheckZeroY true def /PolyEq false def /FXtrans {log} def /FYtrans {Yint sub log} def /RYtrans {10 exch exp Yint add} def /ReduceOrder false def } if \psk@LogTFit\space { /order 1 def /CheckZeroX true def /PolyEq true def /xPr { (log) show (\string\050) show (x) show (\string\051) show } def /FXtrans {log} def /FYtrans {1 mul} def /RYtrans {1 mul} def } if % \psk@LogEFit\space { /order 1 def /CheckZeroX true def /PolyEq true def /xPr { (ln) show (\string\050) show (x) show (\string\051) show } def /FXtrans {ln} def /FYtrans {1 mul} def /RYtrans {1 mul} def } if % \psk@ExpFit\space { /order 1 def /CheckZeroY true def /PolyEq false def /FXtrans {1 mul} def /FYtrans {Yint sub ln} def /RYtrans {Euler exch exp Yint add} def /ReduceOrder false def } if % \psk@GaussFit\space { /order 2 def /CheckZeroY true def /PolyEq false def /FXtrans {1 mul} def /FYtrans {Yint sub ln} def /RYtrans {Euler exch exp Yint add} def /ReduceOrder false def } if % \psk@RecipFit\space { /order 1 def /CheckZeroX true def /PolyEq true def /xPr { (\string\050) show (1) show (\string\244) show (x) show (\string\051) show } def /FXtrans { 1 exch div } def /FYtrans {1 mul} def /RYtrans {1 mul} def } if \psk@CustomFit\space { /PolyFit false def (\psk@FXtrans) length 0 gt { {\psk@FXtrans\space} /FXtrans ED } if (\psk@FYtrans) length 0 gt { {\psk@FYtrans\space} /FYtrans ED } if (\psk@RYtrans) length 0 gt { {\psk@RYtrans\space} /RYtrans ED } if (\psk@CheckZeroX) length 0 gt { {\psk@CheckZeroX\space} /CheckZeroX ED } if (\psk@CheckZeroY) length 0 gt { {\psk@CheckZeroY\space} /CheckZeroY ED } if (\psk@ReduceOrder) length 0 gt { {\psk@ReduceOrder\space} /ReduceOrder ED } if } if % %Initialize Arrays /A order 1 add order 1 add mul array def /B order 1 add array def /ai order 1 add array def %Coordinate Conversion Operators X Y PLU -> Plot Units (pts), X Y RLU -> Real Units (arbitrary) /PLU {exch \pst@number\psxunit\space mul exch \pst@number\psyunit\space mul } def /RLU {exch \pst@number\psxunit\space div exch \pst@number\psyunit\space div } def /PRA { %Inputs: Array [ #rows #cols X Y Scale boolean(print equal sign) boolean(return) ] %"pra" prefix has been added to all variables to minimize chance of naming conflicts % /praSc exch def /praYc exch def /praXc exch def % /praCol exch def /praRow exch def /praARY exch def /praTemp exch def /praARY exch def praTemp length 7 eq {% %praTemp 4 get /praSc exch def %praTemp 3 get /praYc exch def %praTemp 2 get /praXc exch def %praTemp 1 get /praCol exch def %praTemp 0 get /praRow exch def pop /tempfontscale \psk@fontscale praTemp 4 get mul def %/praXc praTemp 2 get \pst@number\psxunit\space mul def /praYc praTemp 3 get \pst@number\psyunit\space mul def /praXmax -999999 def /praYmin 999999 def /prafont {\psk@PSfont\space findfont tempfontscale scalefont setfont} bind def /prasym { \psk@symbolFont\space findfont tempfontscale 2 mul scalefont setfont} bind def prafont praTemp 2 get \pst@number\psxunit\space mul praTemp 3 get \pst@number\psyunit\space mul moveto praARY length praTemp 1 get praTemp 0 get mul eq not { (Row x Col > Vector Dim) show }{% /praXstart praTemp 2 get \pst@number\psxunit\space mul def 0 1 praTemp 1 get 1 sub {/praj ED 0 1 praTemp 0 get 1 sub {/prai ED /prak praj praTemp 1 get prai mul add def praARY prak get \psk@valuewidth string cvs show currentpoint pop dup praXmax gt { /praXmax ED }{ pop } ifelse currentpoint pop praXstart sub neg tempfontscale 1.5 mul neg rmoveto } for currentpoint exch pop dup praYmin lt { /praYmin ED }{ pop } ifelse /praXstart praXmax tempfontscale 0.75 mul add def praXstart praTemp 3 get \pst@number\psyunit\space mul moveto } for %Brackets newpath praTemp 2 get \pst@number\psxunit\space mul tempfontscale 0.125 mul add praTemp 3 get \pst@number\psyunit\space mul tempfontscale add moveto praTemp 2 get \pst@number\psxunit\space mul tempfontscale 0.125 mul sub praTemp 3 get \pst@number\psyunit\space mul tempfontscale add lineto praTemp 2 get \pst@number\psxunit\space mul tempfontscale 0.125 mul sub praYmin tempfontscale add lineto praTemp 2 get \pst@number\psxunit\space mul tempfontscale 0.125 mul add praYmin tempfontscale add lineto stroke newpath praXmax tempfontscale 0.125 mul sub praTemp 3 get \pst@number\psyunit\space mul tempfontscale add moveto praXmax tempfontscale 0.125 mul add praTemp 3 get \pst@number\psyunit\space mul tempfontscale add lineto praXmax tempfontscale 0.125 mul add praYmin tempfontscale add lineto praXmax tempfontscale 0.125 mul sub praYmin tempfontscale add lineto stroke praTemp 5 get { praXmax tempfontscale add praTemp 3 get \pst@number\psyunit\space mul tempfontscale add praYmin add 2 div moveto prasym (=) show currentpoint pop /praXmax ED} if %Left Bracket praTemp 2 get \pst@number\psxunit\space mul praYmin tempfontscale 2 mul sub moveto %Return to bottom left %the following are defined for use with cascade calls to PRA % (Xl,Yb) - Print next matrix below last % (Xr,Yt) - Print next matrix beside last % (Xr,Yb) - Print next matrix beside and below last tx@Dict /praYb known { praYmin tempfontscale 1.5 mul sub \pst@number\psyunit\space div praYb lt { /praYb praYmin tempfontscale 1.5 mul sub \pst@number\psyunit\space div def } if } { /praYb praYmin tempfontscale 1.5 mul sub \pst@number\psyunit\space div def } ifelse /praYb praYmin tempfontscale 1.5 mul sub \pst@number\psyunit\space div def % tx@Dict /praXl known not {/praXl praTemp 2 get def } if % { praTemp 6 get {/praXl praTemp 2 get def} if } ifelse tx@Dict /praYt known not {/praYt praTemp 3 get def } { praTemp 6 get {/praYt praYb def} if } ifelse % tx@Dict /praXr known { praTemp 6 get { /praXr praXl def }{ /praXr praXmax tempfontscale add \pst@number\psxunit\space div def } ifelse }{% /praXr praXmax tempfontscale add \pst@number\psxunit\space div def } ifelse \psk@PSfont\space findfont \psk@fontscale\space scalefont setfont } ifelse }{ \psk@PSfont\space findfont \psk@fontscale\space scalefont setfont 0 0 moveto (PRA Error Not Enough Parameters) show } ifelse } def 0 1 A length 1 sub { A exch 0 put } for 0 1 B length 1 sub { B exch 0 put } for 0 1 ai length 1 sub { ai exch 0 put } for %%%%%%%% %Remove Zero Data (As Required By Transformation Functions) CheckZeroX CheckZeroY or { n { /tmpx false def /tmpy false def CheckZeroY { dup \pst@number\psyunit\space div Yint sub 0 le { /tmpy true def } if } if CheckZeroX { exch dup \pst@number\psxunit\space div 0 le { /tmpx true def } if exch } if tmpx tmpy or { pop pop /n n 1 sub def }{ n 2 mul 2 roll } ifelse } repeat } if %%% %%%% %Determine Max/Min X for Plotting exch dup dup \pst@number\psxunit\space div /xEnd ED \pst@number\psxunit\space div /xStart ED exch /ybar 0 def n { % number of data pairs dup \pst@number\psyunit\space div ybar add /ybar ED 2 copy %copy point \pst@number\psyunit\space div /Yval ED \pst@number\psxunit\space div /Xval ED Xval xStart lt { /xStart Xval def } if % find the lowest xi Xval xEnd gt { /xEnd Xval def } if % find the largest xi n 2 mul 2 roll %roll the original point to the bottom of the stack } repeat ybar n div /ybar ED %Scale X & Y variables and apply Linear Transform if applicable n { \pst@number\psyunit\space div FYtrans n 2 mul 1 roll \pst@number\psxunit\space div FXtrans n 2 mul 1 roll } repeat %START ASSEMBLY OF LEAST SQUARES MATRIX 1 1 order 1 add { /ci exch def 1 1 ci { /cj exch def /epk ci cj add 2 sub def /sum 0 def n { n 2 mul 1 roll dup epk exp sum add /sum ED n 2 mul 1 roll } repeat A ci 1 sub order 1 add mul cj add 1 sub sum put %i row, j col standard A cj 1 sub order 1 add mul ci add 1 sub sum put %j row, i col reversed } for /sum 0 def 1 1 n { /cl exch def 2 copy exch ci 1 sub exp mul sum add /sum ED n 2 mul 2 roll } for B ci 1 sub sum put } for \Pst@Debug\space 1 eq { %A [ order 1 add order 1 add 5 20 0.5 false false ] PRA A [ order 1 add order 1 add \psk@MaPos\space \psk@MaScale\space true false ] PRA B [ order 1 add 1 praXr praYt \psk@MaScale\space false true ] PRA } if %END OF ASSEMBLE MATRIX %REDUCE ORDER FOR Y INTERCEPT FOR MODELS WHICH ARE APPROPRIATE ReduceOrder { /elemOld { exch 1 sub order 1 add mul add 1 sub } def /elemNew { exch 2 sub order mul add 2 sub } def /dummyA order order mul array def 0 1 dummyA length 1 sub { dummyA exch 0 put } for /dummyB order array def 0 1 dummyB length 1 sub { dummyB exch 0 put } for /ai order array def 0 1 ai length 1 sub { ai exch 0 put } for %Translate A to dummyA 2 1 order 1 add {/ci ED 2 1 order 1 add {/cj ED dummyA ci cj elemNew A ci cj elemOld get put } for dummyB ci 2 sub B ci 1 sub get put } for /order order 1 sub def %temporarily reduce system order /A order 1 add order 1 add mul array def 0 1 A length 1 sub {/ci ED A ci dummyA ci get put } for %put dummyA in new A /B order 1 add array def 0 1 B length 1 sub {/ci ED B ci dummyB ci get put } for %put dummyB in new B \Pst@Debug\space 1 eq { A [ order 1 add order 1 add praXr praYt \psk@MaScale\space true false ] PRA B [ order 1 add 1 praXr praYt \psk@MaScale\space false true ] PRA } if \tx@GaussSolve /order order 1 add def %restore system order /m m 1 add def % \Pst@Debug\space 1 eq { %ai [ order 1 add 1 praXr praYt \psk@MaScale\space false false ] PRA } if /dummyai order 1 add array def 0 1 dummyai length 1 sub {/ci ED dummyai ci 0 put } for 0 1 ai length 1 sub {/ci ED dummyai ci 1 add ai ci get put } for /ai dummyai def ai 0 Yint put %ai [ order 2 sub 1 0 0 \psk@MaScale\space true false ] PRA }{ \tx@GaussSolve } ifelse %\Pst@Debug\space 1 eq { %ai [ order 1 add 1 praXr praYt \psk@MaScale\space false false ] PRA } if %%%%%%%%%%%%%%%%%%%%%%%%%% %% END SOLVE CURVE FITS %% %%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%% %% PLOT CURVE FITS %% %%%%%%%%%%%%%%%%%%%%% newpath (\psk@xStart) length 0 gt % special start value? { \psk@xStart\space /xStart ED } if (\psk@xEnd) length 0 gt % special end value? { \psk@xEnd\space /xEnd ED } if /xStep xEnd xStart sub \psk@plotpoints\space div def /Xnow xStart def /Ynow 0 def 1 1 m { /ci exch def ai ci 1 sub get ci 1 sub { Xnow FXtrans mul } repeat Ynow add /Ynow exch def } for Xnow Ynow RYtrans PLU moveto \psk@plotpoints\space cvi { /Xnow Xnow xStep add def /Ynow 0 def 1 1 m { /ci exch def ai ci 1 sub get ci 1 sub { Xnow FXtrans mul } repeat Ynow add /Ynow exch def } for Xnow Ynow RYtrans PLU lineto } repeat % /yminusdot 0 def /yminusbar 0 def n { 2 copy /Ynow ED /Xnow ED /Ydot 0 def 1 1 m { /ci exch def ai ci 1 sub get ci 1 sub { Xnow mul } repeat Ydot add /Ydot exch def } for % Ynow Ydot sub dup mul yminusdot add /yminusdot ED Ynow ybar sub dup mul yminusbar add /yminusbar ED n 2 mul 2 roll } repeat 1 yminusdot yminusbar div sub /rsquared ED %%%%%%%%%%%%%%%%%%%%%%%%% %% END PLOT CURVE FITS %% %%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%% %% START PRINT EQUATIONS %% %%%%%%%%%%%%%%%%%%%%%%%%%%% \ifPst@ShowEq % Setup Cursor (\psk@EqPos) length 0 gt {\psk@EqPos\space PLU } {2 -5 PLU } ifelse moveto % Power Fit Equation % \psk@PowerFit\space { 10 ai 0 get exp \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /tmpA exch def ai 1 get \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /tmpB exch def (y = ) show tmpA \psk@valuewidth string cvs show (x) show 0 \psk@fontscale 4 div rmoveto mfontexp tmpB \psk@valuewidth string cvs show mfont 0 \psk@fontscale 4 div neg rmoveto setYint Yint 0 ne and { Yint 0 gt { (+) show } if Yint \psk@valuewidth string cvs show } if } if % Exponential Fit Equation % \psk@ExpFit\space { Euler ai 0 get exp \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /tmpA exch def Euler ai 1 get exp \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /tmpB exch def (y = ) show (\string\050) show tmpA \psk@valuewidth string cvs show (\string\051) show %(e) show (\string\264) show %0 \psk@fontscale 4 div rmoveto %mfontexp tmpB \psk@valuewidth string cvs show (\string\050) show tmpB \psk@valuewidth string cvs show (\string\051) show 0 \psk@fontscale 2 div rmoveto mfontexp (x) show mfont 0 \psk@fontscale 2 div neg rmoveto setYint Yint 0 ne and { Yint 0 gt { (+) show } if Yint \psk@valuewidth string cvs show } if} if % Gaussian Fit Equation % \psk@GaussFit\space { /tmpA Euler ai 0 get ai 1 get dup mul 4 div ai 2 get div sub exp def /tmpB ai 1 get 2 div ai 2 get div def /tmpC ai 2 get def tmpA \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /tmpA exch def tmpB \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /tmpB exch def tmpC \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /tmpC exch def (y = ) show tmpA 30 string cvs show (e) show 0 \psk@fontscale 4 div rmoveto mfontexp tmpC \psk@valuewidth string cvs show (\string\(x) show tmpB 0 ne { tmpB 0 ge { (+) show } if tmpB \psk@valuewidth string cvs show } if (\string\)) show 0 \psk@fontscale 4 div rmoveto (2) show mfont 0 \psk@fontscale 2 div neg rmoveto setYint Yint 0 ne and { Yint 0 gt { (+) show } if Yint \psk@valuewidth string cvs show } if } if %%%%%%%%%%%%%%% % Round Coeff % %%%%%%%%%%%%%%% 0 1 ai length 1 sub {/ci exch def ai ci get \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if ai ci 3 -1 roll put } for %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Plot Polynomial (including Linear) Fit % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% PolyEq { (y =) show 1 1 m { /ci exch def ai ci 1 sub get 30 string cvs show ci 1 sub 0 gt { xPr } if ci 1 sub 1 gt {% 0 \psk@fontscale 4 div rmoveto ci 1 sub 30 string cvs show 0 \psk@fontscale 4 div neg rmoveto } if ci m eq not { ai ci get 0 lt {}{ (+) show } ifelse } if } for } if %Print Rsquared rsquared \psk@decimals\space -1 gt { 10 \psk@decimals exp dup 3 1 roll mul round exch div } if \psk@decimals\space 0 eq { cvi } if /rsquared exch def \psk@EqPos\space PLU moveto 0 \psk@fontscale 1.5 mul neg rmoveto (R) show 0 \psk@fontscale 2 div rmoveto mfontexp (2) show mfont 0 \psk@fontscale 2 div neg rmoveto ( = ) show rsquared \psk@valuewidth string cvs show \fi } % \def\tx@GaussSolve{GaussSolve } \pst@def{GaussSolve}<{% %Required Variables %order,A,B,ai /er 1 def % /tol 0.00001 def % /m order 1 add def %number of equations /Aout m m mul array def %This is for printing the scaled identity matrix which A becomes 0 1 Aout length 1 sub { Aout exch 0 put } for /elem { exch 1 sub m mul add 1 sub } def %Array i j %DETERMINE SCALE FACTORS /sc m array def %row scale factors 0 1 sc length 1 sub { sc exch 0 put } for 1 1 m { /ci exch def sc ci 1 sub A ci 1 elem get dup mul sqrt put 2 1 m { /cj exch def A ci cj elem get dup mul sqrt sc ci 1 sub get gt {% sc ci 1 sub A ci cj elem get dup mul sqrt put } if } for % 1 ci mul 1 sc ci 1 sub get 30 string cvs show } for %END OF SCALE FACTORS %CALL ELIMINATE order 1 ge {% 1 1 m 1 sub { /ck exch def %Call Pivot /p ck def A ck ck elem get sc ck 1 sub get div dup mul sqrt /big exch def ck 1 add 1 m {/cii exch def A cii ck elem get sc cii 1 sub get div dup mul sqrt /dummy exch def dummy big gt { /big dummy def /p cii def } if } for p ck eq not { ck 1 m {/cjj exch def A p cjj elem get /dummy exch def A p cjj elem A ck cjj elem get put A ck cjj elem dummy put } for B p 1 sub get /dummy exch def B p 1 sub B ck 1 sub get put B ck 1 sub dummy put sc p 1 sub get /dummy exch def sc p 1 sub sc ck 1 sub get put sc ck 1 sub dummy put } if %End Pivot %Back to Eliminate A ck ck elem get sc ck 1 sub get div dup mul sqrt tol lt { /er -1 def } if ck 1 add 1 m { /ci exch def A ci ck elem get A ck ck elem get div /factor exch def %ck 1 add 1 m { /cj exch def 1 1 m { /cj exch def %Line Changed from above, requires more calcs but allows debugging of matrix ops. A ci cj elem %index A ci cj elem get A ck cj elem get factor mul sub %value put } for B ci 1 sub %index B ci 1 sub get B ck 1 sub get factor mul sub %value put } for } for \Pst@Debug\space 1 eq { A [ order 1 add order 1 add praXr praYt \psk@MaScale\space true false ] PRA B [ order 1 add 1 praXr praYt \psk@MaScale\space false true ] PRA } if } if %Start Back Substitution %Add Error check for er -1 ai m 1 sub B m 1 sub get A m m elem get div put %%%%Added for debug Aout m m elem A m m elem get dup div put %%%% order 1 ge {% m 1 sub -1 1 { /ci exch def /sum 0 def ci 1 add 1 m { /cj exch def A ci cj elem get ai cj 1 sub get mul sum add /sum exch def } for ai ci 1 sub B ci 1 sub get sum sub A ci ci elem get div put %%%%Added for debug Aout ci ci elem A ci ci elem get dup div put %%%% } for A ck ck elem get sc ck 1 sub get div dup mul sqrt tol lt { /er -1 def } if } if \Pst@Debug\space 1 eq { Aout [ order 1 add order 1 add praXr praYt \psk@MaScale\space true false ] PRA ai [ order 1 add 1 praXr praYt \psk@MaScale\space false false ] PRA } if}> % \catcode`\@=\PstAtCode\relax % %% END: pst-fit.tex \endinput