%! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % POSTSCRIPT LINEAR EQUATION SOLVER USING DETERMINANTS % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Version 1.3 September 15, 1996. % Copyright c 1996 by Don Lancaster and Synergetics 3860 West First % Street, Box 809, Thatcher AZ, 85552. don@tinaja.com Full technical % support on GURU'S LAIR at www.tinaja.com % All commercial rights and all electronic media rights are **FULLY** % reserved. Reposting expressly forbidden. Consulting services available. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % This is a continuing example in our series of using PostScript as % a general purpose computing language. Two way host recordable com % is **REQUIRED** for portions of these utilities. % View the file initially with an editor or word processor. Then comment % or uncomment the various demos. Finally, extract or customize the % actual PostScript code as stand-alone utilities. These utilities and % demos may be viewed with GhostScript, or routed to any PostScript device. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Digital filters are but one of many places where high variable % linear algebraic equations need to be solved. The method of % determinants gives us an orderly way to solve, say, ten equations % in ten unknowns in two seconds. % In this novel PS code, a "modified" determinant method is used for faster % solutions. All but one element of each determinant row is forced % to zero. Then a determinant of any given order is solved for only % one variable. The result is used in turn to solve the next lower order. % More details on determinants themselves and speedup techniques % appear in MUSE106.PDF on www.tinaja.com Also in the December? % 1996 issue of Electronics Now, and in the TECH MUSINGS reprints % from SYNERGETICS. % This file first shows how to solve determinants, then linear % equations, then detailed linear phase filter response problems. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % COMPUTE A BRUTE FORCE 4X4 DETERMINANT % ===================================== % /det4 accepts a b c d e f g h i j k l m n o p on the stack and % computes the resultant 4 x 4 determinant def's trades vm for speed. /det4 { % save /det4snap exch def /d4p exch def /d4o exch def /d4n exch def /d4m exch def /d4l exch def /d4k exch def /d4j exch def /d4i exch def /d4h exch def /d4g exch def /d4f exch def /d4e exch def /d4d exch def /d4c exch def /d4b exch def /d4a exch def d4a d4f d4k d4p mul d4l d4o mul sub mul d4g d4j d4p mul d4l d4n mul sub mul sub d4h d4j d4o mul d4k d4n mul sub mul add mul d4b d4e d4k d4p mul d4l d4o mul sub mul d4g d4i d4p mul d4l d4m mul sub mul sub d4h d4i d4o mul d4k d4m mul sub mul add mul sub d4c d4e d4j d4p mul d4l d4n mul sub mul d4f d4i d4p mul d4l d4m mul sub mul sub d4h d4i d4n mul d4m d4j mul sub mul add mul add d4d d4e d4j d4o mul d4n d4k mul sub mul d4f d4i d4o mul d4k d4m mul sub mul sub d4g d4i d4n mul d4j d4m mul sub mul add mul sub } bind def %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % COMPUTE A FAST 5X5 DETERMINANT % ============================== % /det5 accepts a b c d e f g h .... w x y on the stack and % computes the resultant 5 x 5 determinant. Four zeros are forced. /det5 {save /det5snap exch def /d5y exch def /d5x exch def /d5w exch def /d5v exch def /d5u exch def /d5t exch def /d5s exch def /d5r exch def /d5q exch def /d5p exch def /d5o exch def /d5n exch def /d5m exch def /d5l exch def /d5k exch def /d5j exch def /d5i exch def /d5h exch def /d5g exch def /d5f exch def /d5e exch def /d5d exch def /d5c exch def /d5b exch def /d5a exch def % a b c d e ----> a b c d (ma + e) force (ma + e) = 0 m = -e/a % f g h i j % k l m n o % p q r s t % u v w x y d5a 0 ne { /mb d5b d5a div neg def /d5b 0 store /d5g d5f mb mul d5g add store /d5l d5k mb mul d5l add store /d5q d5p mb mul d5q add store /d5v d5u mb mul d5v add store } if d5a 0 ne { /mc d5c d5a div neg def /d5c 0 store /d5h d5f mc mul d5h add store /d5m d5k mc mul d5m add store /d5r d5p mc mul d5r add store /d5w d5u mc mul d5w add store } if d5a 0 ne { /md d5d d5a div neg def /d5d 0 store /d5i d5f md mul d5i add store /d5n d5k md mul d5n add store /d5s d5p md mul d5s add store /d5x d5u md mul d5x add store } if d5a 0 ne { /me d5e d5a div neg def /d5e 0 store /d5j d5f me mul d5j add store /d5o d5k me mul d5o add store /d5t d5p me mul d5t add store /d5y d5u me mul d5y add store } if d5a d5g d5h d5i d5j d5l d5m d5n d5o d5q d5r d5s d5t d5v d5w d5x d5y det4 mul d5b 0 ne { d5b d5f d5h d5i d5j d5k d5m d5n d5o d5p d5r d5s d5t d5u d5w d5x d5y det4 mul sub } if d5c 0 ne { d5c d5f d5g d5i d5j d5k d5l d5n d5o d5p d5q d5s d5t d5u d5v d5x d5y det4 mul add } if d5d 0 ne { d5d d5f d5g d5h d5j d5k d5l d5m d5o d5p d5q d5r d5t d5u d5v d5w d5y det4 mul sub } if d5e 0 ne { d5e d5f d5g d5h d5i d5k d5l d5m d5n d5p d5q d5r d5s d5u d5v d5w d5x det4 mul add } if det5snap restore } def % For instance 1 3 2 4 5 6 7 8 9 1 3 2 5 4 3 2 % 1 2 4 3 1 3 4 2 1 det5 is 158 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % COMPUTE A FAST 6X6 DETERMINANT % ============================== % /det6 accepts a b c d e f g h .... w x y ... ii jj on the stack and % computes the resultant 6 x 6 determinant. Five zeros are forced /det6 {save /det6snap exch def /d6jj exch def /d6ii exch def /d6hh exch def /d6gg exch def /d6ff exch def /d6ee exch def /d6dd exch def /d6cc exch def /d6bb exch def /d6aa exch def /d6z exch def /d6y exch def /d6x exch def /d6w exch def /d6v exch def /d6u exch def /d6t exch def /d6s exch def /d6r exch def /d6q exch def /d6p exch def /d6o exch def /d6n exch def /d6m exch def /d6l exch def /d6k exch def /d6j exch def /d6i exch def /d6h exch def /d6g exch def /d6f exch def /d6e exch def /d6d exch def /d6c exch def /d6b exch def /d6a exch def % a b c d e f ----> a b c d (ma + e) force (ma + e) = 0 m = -e/a % g h i j k l % m n o p q r % s t u v w x % y z A B C D % E F G H I J d6a 0 ne { /mf d6f d6a div neg def /d6f 0 store /d6l d6g mf mul d6l add store /d6r d6m mf mul d6r add store /d6x d6s mf mul d6x add store /d6dd d6y mf mul d6dd add store /d6jj d6ee mf mul d6jj add store } if d6a 0 ne { /me d6e d6a div neg def /d6e 0 store /d6k d6g me mul d6k add store /d6q d6m me mul d6q add store /d6w d6s me mul d6w add store /d6cc d6y me mul d6cc add store /d6ii d6ee me mul d6ii add store } if d6a 0 ne { /md d6d d6a div neg def /d6d 0 store /d6j d6g md mul d6j add store /d6p d6m md mul d6p add store /d6v d6s md mul d6v add store /d6bb d6y md mul d6bb add store /d6hh d6ee md mul d6hh add store } if d6a 0 ne { /mc d6c d6a div neg def /d6c 0 store /d6i d6g mc mul d6i add store /d6o d6m mc mul d6o add store /d6u d6s mc mul d6u add store /d6aa d6y mc mul d6aa add store /d6gg d6ee mc mul d6gg add store } if d6a 0 ne { /mb d6b d6a div neg def /d6b 0 store /d6h d6g mb mul d6h add store /d6n d6m mb mul d6n add store /d6t d6s mb mul d6t add store /d6z d6y mb mul d6z add store /d6ff d6ee mb mul d6ff add store } if d6a d6h d6i d6j d6k d6l d6n d6o d6p d6q d6r d6t d6u d6v d6w d6x d6z d6aa d6bb d6cc d6dd d6ff d6gg d6hh d6ii d6jj det5 mul d6b 0 ne { d6b d6g d6i d6j d6k d6l d6m d6o d6p d6q d6r d6s d6u d6v d6w d6x d6y d6aa d6bb d6cc d6dd d6ee d6gg d6hh d6ii d6jj det5 mul sub} if d6c 0 ne { d6c d6g d6h d6j d6k d6l d6m d6n d6p d6q d6r d6s d6t d6v d6w d6x d6y d6z d6bb d6cc d6dd d6ee d6ff d6hh d6ii d6jj det5 mul add} if d6d 0 ne { d6d d6g d6h d6i d6k d6l d6m d6n d6o d6q d6r d6s d6t d6u d6w d6x d6y d6z d6aa d6cc d6dd d6ee d6ff d6gg d6ii d6jj det5 mul sub} if d6e 0 ne { d6e d6g d6h d6i d6j d6l d6m d6n d6o d6p d6r d6s d6t d6u d6v d6x d6y d6z d6aa d6bb d6dd d6ee d6ff d6gg d6hh d6jj det5 mul add} if d6f 0 ne { d6f d6g d6h d6i d6j d6k d6m d6n d6o d6p d6q d6s d6t d6u d6v d6w d6y d6z d6aa d6bb d6cc d6ee d6ff d6gg d6hh d6ii det5 mul sub } if det6snap restore } bind def % For instance % 1 3 2 4 5 6 % 7 8 9 1 3 2 % 5 4 3 2 1 2 % 4 3 1 3 4 2 % 1 2 5 3 7 2 % 5 3 3 1 2 3 det6 is a -2976 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % COMPUTE A 7X7 FAST DETERMINANT % ============================== % /det7 accepts a b c d e f g h .... w x y ... vv ww on the stack and % computes the resultant 7 x 7 determinant. Six zeros are forced. /det7 {save /det7snap exch def /d7ww exch def /d7vv exch def /d7uu exch def /d7tt exch def /d7ss exch def /d7rr exch def /d7qq exch def /d7pp exch def /d7oo exch def /d7nn exch def /d7mm exch def /d7ll exch def /d7kk exch def /d7jj exch def /d7ii exch def /d7hh exch def /d7gg exch def /d7ff exch def /d7ee exch def /d7dd exch def /d7cc exch def /d7bb exch def /d7aa exch def /d7z exch def /d7y exch def /d7x exch def /d7w exch def /d7v exch def /d7u exch def /d7t exch def /d7s exch def /d7r exch def /d7q exch def /d7p exch def /d7o exch def /d7n exch def /d7m exch def /d7l exch def /d7k exch def /d7j exch def /d7i exch def /d7h exch def /d7g exch def /d7f exch def /d7e exch def /d7d exch def /d7c exch def /d7b exch def /d7a exch def % a b c d e f g ----> a b c d (ma + e) force (ma + e) = 0 m = -e/a % h i j k l m n % o p q r s t u % v w x y z A B % C D E F G H I % J K L M N O P % Q R S T U V W d7a 0 ne { /mb d7b d7a div neg def /d7b 0 store /d7i d7h mb mul d7i add store /d7p d7o mb mul d7p add store /d7w d7v mb mul d7w add store /d7dd d7cc mb mul d7dd add store /d7kk d7jj mb mul d7kk add store /d7rr d7qq mb mul d7rr add store } if d7a 0 ne { /mc d7c d7a div neg def /d7c 0 store /d7j d7h mc mul d7j add store /d7q d7o mc mul d7q add store /d7x d7v mc mul d7x add store /d7ee d7cc mc mul d7ee add store /d7ll d7jj mc mul d7ll add store /d7ss d7qq mc mul d7ss add store } if d7a 0 ne { /md d7d d7a div neg def /d7d 0 store /d7k d7h md mul d7k add store /d7r d7o md mul d7r add store /d7y d7v md mul d7y add store /d7ff d7cc md mul d7ff add store /d7mm d7jj md mul d7mm add store /d7tt d7qq md mul d7tt add store } if d7a 0 ne { /me d7e d7a div neg def /d7e 0 store /d7l d7h me mul d7l add store /d7s d7o me mul d7s add store /d7z d7v me mul d7z add store /d7gg d7cc me mul d7gg add store /d7nn d7jj me mul d7nn add store /d7uu d7qq me mul d7uu add store } if d7a 0 ne { /mf d7f d7a div neg def /d7f 0 store /d7m d7h mf mul d7m add store /d7t d7o mf mul d7t add store /d7aa d7v mf mul d7aa add store /d7hh d7cc mf mul d7hh add store /d7oo d7jj mf mul d7oo add store /d7vv d7qq mf mul d7vv add store } if d7a 0 ne { /mg d7g d7a div neg def /d7g 0 store /d7n d7h mg mul d7n add store /d7u d7o mg mul d7u add store /d7bb d7v mg mul d7bb add store /d7ii d7cc mg mul d7ii add store /d7pp d7jj mg mul d7pp add store /d7ww d7qq mg mul d7ww add store } if d7a 0 ne { d7a d7i d7j d7k d7l d7m d7n d7p d7q d7r d7s d7t d7u d7w d7x d7y d7z d7aa d7bb d7dd d7ee d7ff d7gg d7hh d7ii d7kk d7ll d7mm d7nn d7oo d7pp d7rr d7ss d7tt d7uu d7vv d7ww det6 mul}{0} ifelse d7b 0 ne { d7b d7h d7j d7k d7l d7m d7n d7o d7q d7r d7s d7t d7u d7v d7x d7y d7z d7aa d7bb d7cc d7ee d7ff d7gg d7hh d7ii d7jj d7ll d7mm d7nn d7oo d7pp d7qq d7ss d7tt d7uu d7vv d7ww det6 mul sub} if d7c 0 ne { d7c d7h d7i d7k d7l d7m d7n d7o d7p d7r d7s d7t d7u d7v d7w d7y d7z d7aa d7bb d7cc d7dd d7ff d7gg d7hh d7ii d7jj d7kk d7mm d7nn d7oo d7pp d7qq d7rr d7tt d7uu d7vv d7ww det6 mul add} if d7d 0 ne { d7d d7h d7i d7j d7l d7m d7n d7o d7p d7q d7s d7t d7u d7v d7w d7x d7z d7aa d7bb d7cc d7dd d7ee d7gg d7hh d7ii d7jj d7kk d7ll d7nn d7oo d7pp d7qq d7rr d7ss d7uu d7vv d7ww det6 mul sub} if d7e 0 ne { d7e d7h d7i d7j d7k d7m d7n d7o d7p d7q d7r d7t d7u d7v d7w d7x d7y d7aa d7bb d7cc d7dd d7ee d7ff d7hh d7ii d7jj d7kk d7ll d7mm d7oo d7pp d7qq d7rr d7ss d7tt d7vv d7ww det6 mul add } if d7f 0 ne { d7f d7h d7i d7j d7k d7l d7n d7o d7p d7q d7r d7s d7u d7v d7w d7x d7y d7z d7bb d7cc d7dd d7ee d7ff d7gg d7ii d7jj d7kk d7ll d7mm d7nn d7pp d7qq d7rr d7ss d7tt d7uu d7ww det6 mul sub} if d7g 0 ne { d7g d7h d7i d7j d7k d7l d7m d7o d7p d7q d7r d7s d7t d7v d7w d7x d7y d7z d7aa d7cc d7dd d7ee d7ff d7gg d7hh d7jj d7kk d7ll d7mm d7nn d7oo d7qq d7rr d7ss d7tt d7uu d7vv det6 mul add } if det7snap restore } bind def %%%%%%%%%%%%%%%% % COMPUTE A 8X8 FAST DETERMINANT % ============================== % /det8 accepts a b .... y z aa ab ... bl bm on the stack and % computes the resultant 8 x 8 determinant. Seven zeros are forced /det8 {save /det8snap exch def /d8bl exch def /d8bk exch def /d8bj exch def /d8bi exch def /d8bh exch def /d8bg exch def /d8bf exch def /d8be exch def /d8bd exch def /d8bc exch def /d8bb exch def /d8ba exch def /d8az exch def /d8ay exch def /d8ax exch def /d8aw exch def /d8av exch def /d8au exch def /d8at exch def /d8as exch def /d8ar exch def /d8aq exch def /d8ap exch def /d8ao exch def /d8an exch def /d8am exch def /d8al exch def /d8ak exch def /d8aj exch def /d8ai exch def /d8ah exch def /d8ag exch def /d8af exch def /d8ae exch def /d8ad exch def /d8ac exch def /d8ab exch def /d8aa exch def /d8z exch def /d8y exch def /d8x exch def /d8w exch def /d8v exch def /d8u exch def /d8t exch def /d8s exch def /d8r exch def /d8q exch def /d8p exch def /d8o exch def /d8n exch def /d8m exch def /d8l exch def /d8k exch def /d8j exch def /d8i exch def /d8h exch def /d8g exch def /d8f exch def /d8e exch def /d8d exch def /d8c exch def /d8b exch def /d8a exch def % a b c d e f g h ----> a b c d (ma + e) force (ma + e) = 0 m = -e/a % i j k l m n o p % q r s t u v w x % y z A B C D E F % G H I J K L M N % O P Q R S T U V % W X Y Z a b c d % e f g h i j k l d8a 0 ne { /mh d8h d8a div neg def /d8h 0 store /d8p d8i mh mul d8p add store /d8x d8q mh mul d8x add store /d8af d8y mh mul d8af add store /d8an d8ag mh mul d8an add store /d8av d8ao mh mul d8av add store /d8bd d8aw mh mul d8bd add store /d8bl d8be mh mul d8bl add store } if d8a 0 ne { /mg d8g d8a div neg def /d8g 0 store /d8o d8i mg mul d8o add store /d8w d8q mg mul d8w add store /d8ae d8y mg mul d8ae add store /d8am d8ag mg mul d8am add store /d8au d8ao mg mul d8au add store /d8bc d8aw mg mul d8bc add store /d8bk d8be mg mul d8bk add store } if d8a 0 ne { /mf d8f d8a div neg def /d8f 0 store /d8n d8i mf mul d8n add store /d8v d8q mf mul d8v add store /d8ad d8y mf mul d8ad add store /d8al d8ag mf mul d8al add store /d8at d8ao mf mul d8at add store /d8bb d8aw mf mul d8bb add store /d8bj d8be mf mul d8bj add store } if d8a 0 ne { /me d8e d8a div neg def /d8e 0 store /d8m d8i me mul d8m add store /d8u d8q me mul d8u add store /d8ac d8y me mul d8ac add store /d8ak d8ag me mul d8ak add store /d8as d8ao me mul d8as add store /d8ba d8aw me mul d8ba add store /d8bi d8be me mul d8bi add store } if d8a 0 ne { /md d8d d8a div neg def /d8d 0 store /d8l d8i md mul d8l add store /d8t d8q md mul d8t add store /d8ab d8y md mul d8ab add store /d8aj d8ag md mul d8aj add store /d8ar d8ao md mul d8ar add store /d8az d8aw md mul d8az add store /d8bh d8be md mul d8bh add store } if d8a 0 ne { /mc d8c d8a div neg def /d8c 0 store /d8k d8i mc mul d8k add store /d8s d8q mc mul d8s add store /d8aa d8y mc mul d8aa add store /d8ai d8ag mc mul d8ai add store /d8aq d8ao mc mul d8aq add store /d8ay d8aw mc mul d8ay add store /d8bg d8be mc mul d8bg add store } if d8a 0 ne { /mb d8b d8a div neg def /d8b 0 store /d8j d8i mb mul d8j add store /d8r d8q mb mul d8r add store /d8z d8y mb mul d8z add store /d8ah d8ag mb mul d8ah add store /d8ap d8ao mb mul d8ap add store /d8ax d8aw mb mul d8ax add store /d8bf d8be mb mul d8bf add store } if d8a 0 ne { d8a d8j d8k d8l d8m d8n d8o d8p d8r d8s d8t d8u d8v d8w d8x d8z d8aa d8ab d8ac d8ad d8ae d8af d8ah d8ai d8aj d8ak d8al d8am d8an d8ap d8aq d8ar d8as d8at d8au d8av d8ax d8ay d8az d8ba d8bb d8bc d8bd d8bf d8bg d8bh d8bi d8bj d8bk d8bl det7 mul % dup report } {0} ifelse d8b 0 ne { d8b d8i d8k d8l d8m d8n d8o d8p d8q d8s d8t d8u d8v d8w d8x d8y d8aa d8ab d8ac d8ad d8ae d8af d8ag d8ai d8aj d8ak d8al d8am d8an d8ao d8aq d8ar d8as d8at d8au d8av d8aw d8ay d8az d8ba d8bb d8bc d8bd d8be d8bg d8bh d8bi d8bj d8bk d8bl det7 mul sub} if d8c 0 ne { d8c d8i d8j d8l d8m d8n d8o d8p d8q d8r d8t d8u d8v d8w d8x d8y d8z d8ab d8ac d8ad d8ae d8af d8ag d8ah d8aj d8ak d8al d8am d8an d8ao d8ap d8ar d8as d8at d8au d8av d8aw d8ax d8az d8ba d8bb d8bc d8bd d8be d8bf d8bh d8bi d8bj d8bk d8bl det7 mul add} if d8d 0 ne { d8d d8i d8j d8k d8m d8n d8o d8p d8q d8r d8s d8u d8v d8w d8x d8y d8z d8aa d8ac d8ad d8ae d8af d8ag d8ah d8ai d8ak d8al d8am d8an d8ao d8ap d8aq d8as d8at d8au d8av d8aw d8ax d8ay d8ba d8bb d8bc d8bd d8be d8bf d8bg d8bi d8bj d8bk d8bl det7 mul sub} if d8e 0 ne { d8e d8i d8j d8k d8l d8n d8o d8p d8q d8r d8s d8t d8v d8w d8x d8y d8z d8aa d8ab d8ad d8ae d8af d8ag d8ah d8ai d8aj d8al d8am d8an d8ao d8ap d8aq d8ar d8at d8au d8av d8aw d8ax d8ay d8az d8bb d8bc d8bd d8be d8bf d8bg d8bh d8bj d8bk d8bl det7 mul add} if d8f 0 ne { d8f d8i d8j d8k d8l d8m d8o d8p d8q d8r d8s d8t d8u d8w d8x d8y d8z d8aa d8ab d8ac d8ae d8af d8ag d8ah d8ai d8aj d8ak d8am d8an d8ao d8ap d8aq d8ar d8as d8au d8av d8aw d8ax d8ay d8az d8ba d8bc d8bd d8be d8bf d8bg d8bh d8bi d8bk d8bl det7 mul sub} if d8g 0 ne { d8g d8i d8j d8k d8l d8m d8n d8p d8q d8r d8s d8t d8u d8v d8x d8y d8z d8aa d8ab d8ac d8ad d8af d8ag d8ah d8ai d8aj d8ak d8al d8an d8ao d8ap d8aq d8ar d8as d8at d8av d8aw d8ax d8ay d8az d8ba d8bb d8bd d8be d8bf d8bg d8bh d8bi d8bj d8bl det7 mul add} if d8h 0 ne { d8h d8i d8j d8k d8l d8m d8n d8o d8q d8r d8s d8t d8u d8v d8w d8y d8z d8aa d8ab d8ac d8ad d8ae d8ag d8ah d8ai d8aj d8ak d8al d8am d8ao d8ap d8aq d8ar d8as d8at d8au d8aw d8ax d8ay d8az d8ba d8bb d8bc d8be d8bf d8bg d8bh d8bi d8bj d8bk det7 mul sub} if det8snap restore } bind def %%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % COMPUTE A 9X9 FAST DETERMINANT % ============================== % /det9 accepts a b .... y z aa ab ... cb cc on the stack and % computes the resultant 9 x 9 determinant. Eight zeros are forced. /det9 {save /det9snap exch def /d9cc exch def /d9cb exch def /d9ca exch def /d9bz exch def /d9by exch def /d9bx exch def /d9bw exch def /d9bv exch def /d9bu exch def /d9bt exch def /d9bs exch def /d9br exch def /d9bq exch def /d9bp exch def /d9bo exch def /d9bn exch def /d9bm exch def /d9bl exch def /d9bk exch def /d9bj exch def /d9bi exch def /d9bh exch def /d9bg exch def /d9bf exch def /d9be exch def /d9bd exch def /d9bc exch def /d9bb exch def /d9ba exch def /d9az exch def /d9ay exch def /d9ax exch def /d9aw exch def /d9av exch def /d9au exch def /d9at exch def /d9as exch def /d9ar exch def /d9aq exch def /d9ap exch def /d9ao exch def /d9an exch def /d9am exch def /d9al exch def /d9ak exch def /d9aj exch def /d9ai exch def /d9ah exch def /d9ag exch def /d9af exch def /d9ae exch def /d9ad exch def /d9ac exch def /d9ab exch def /d9aa exch def /d9z exch def /d9y exch def /d9x exch def /d9w exch def /d9v exch def /d9u exch def /d9t exch def /d9s exch def /d9r exch def /d9q exch def /d9p exch def /d9o exch def /d9n exch def /d9m exch def /d9l exch def /d9k exch def /d9j exch def /d9i exch def /d9h exch def /d9g exch def /d9f exch def /d9e exch def /d9d exch def /d9c exch def /d9b exch def /d9a exch def % a b c d e f g h i % j k l m n o p q r % s t u v w x y z A % B C D E F G H I J % K L M N O P Q R S % T U V W X Y Z a b % c d e f g h i j k % l m n o p q r s t % u v w x y z A B C d9a 0 ne { /mb d9b d9a div neg def /d9b 0 store /d9k d9j mb mul d9k add store /d9t d9s mb mul d9t add store /d9ac d9ab mb mul d9ac add store /d9al d9ak mb mul d9al add store /d9au d9at mb mul d9au add store /d9bd d9bc mb mul d9bd add store /d9bm d9bl mb mul d9bm add store /d9bv d9bu mb mul d9bv add store } if d9a 0 ne { /mc d9c d9a div neg def /d9c 0 store /d9l d9j mc mul d9l add store /d9u d9s mc mul d9u add store /d9ad d9ab mc mul d9ad add store /d9am d9ak mc mul d9am add store /d9av d9at mc mul d9av add store /d9be d9bc mc mul d9be add store /d9bn d9bl mc mul d9bn add store /d9bw d9bu mc mul d9bw add store } if d9a 0 ne { /md d9d d9a div neg def /d9d 0 store /d9m d9j md mul d9m add store /d9v d9s md mul d9v add store /d9ae d9ab md mul d9ae add store /d9an d9ak md mul d9an add store /d9aw d9at md mul d9aw add store /d9bf d9bc md mul d9bf add store /d9bo d9bl md mul d9bo add store /d9bx d9bu md mul d9bx add store } if d9a 0 ne { /me d9e d9a div neg def /d9e 0 store /d9n d9j me mul d9n add store /d9w d9s me mul d9w add store /d9af d9ab me mul d9af add store /d9ao d9ak me mul d9ao add store /d9ax d9at me mul d9ax add store /d9bg d9bc me mul d9bg add store /d9bp d9bl me mul d9bp add store /d9by d9bu me mul d9by add store } if d9a 0 ne { /mf d9f d9a div neg def /d9f 0 store /d9o d9j mf mul d9o add store /d9x d9s mf mul d9x add store /d9ag d9ab mf mul d9ag add store /d9ap d9ak mf mul d9ap add store /d9ay d9at mf mul d9ay add store /d9bh d9bc mf mul d9bh add store /d9bq d9bl mf mul d9bq add store /d9bz d9bu mf mul d9bz add store } if d9a 0 ne { /mg d9g d9a div neg def /d9g 0 store /d9p d9j mg mul d9p add store /d9y d9s mg mul d9y add store /d9ah d9ab mg mul d9ah add store /d9aq d9ak mg mul d9aq add store /d9az d9at mg mul d9az add store /d9bi d9bc mg mul d9bi add store /d9br d9bl mg mul d9br add store /d9ca d9bu mg mul d9ca add store } if d9a 0 ne { /mh d9h d9a div neg def /d9h 0 store /d9q d9j mh mul d9q add store /d9z d9s mh mul d9z add store /d9ai d9ab mh mul d9ai add store /d9ar d9ak mh mul d9ar add store /d9ba d9at mh mul d9ba add store /d9bj d9bc mh mul d9bj add store /d9bs d9bl mh mul d9bs add store /d9cb d9bu mh mul d9cb add store } if d9a 0 ne { /mi d9i d9a div neg def /d9i 0 store /d9r d9j mi mul d9r add store /d9aa d9s mi mul d9aa add store /d9aj d9ab mi mul d9aj add store /d9as d9ak mi mul d9as add store /d9bb d9at mi mul d9bb add store /d9bk d9bc mi mul d9bk add store /d9bt d9bl mi mul d9bt add store /d9cc d9bu mi mul d9cc add store } if d9a 0 ne { d9a d9k d9l d9m d9n d9o d9p d9q d9r d9t d9u d9v d9w d9x d9y d9z d9aa d9ac d9ad d9ae d9af d9ag d9ah d9ai d9aj d9al d9am d9an d9ao d9ap d9aq d9ar d9as d9au d9av d9aw d9ax d9ay d9az d9ba d9bb d9bd d9be d9bf d9bg d9bh d9bi d9bj d9bk d9bm d9bn d9bo d9bp d9bq d9br d9bs d9bt d9bv d9bw d9bx d9by d9bz d9ca d9cb d9cc det8 mul }{0} ifelse d9b 0 ne { d9b d9j d9l d9m d9n d9o d9p d9q d9r d9s d9u d9v d9w d9x d9y d9z d9aa d9ab d9ad d9ae d9af d9ag d9ah d9ai d9aj d9ak d9am d9an d9ao d9ap d9aq d9ar d9as d9at d9av d9aw d9ax d9ay d9az d9ba d9bb d9bc d9be d9bf d9bg d9bh d9bi d9bj d9bk d9bl d9bn d9bo d9bp d9bq d9br d9bs d9bt d9bu d9bw d9bx d9by d9bz d9ca d9cb d9cc det8 mul sub } if d9c 0 ne { d9c d9j d9k d9m d9n d9o d9p d9q d9r d9s d9t d9v d9w d9x d9y d9z d9aa d9ab d9ac d9ae d9af d9ag d9ah d9ai d9aj d9ak d9al d9an d9ao d9ap d9aq d9ar d9as d9at d9au d9aw d9ax d9ay d9az d9ba d9bb d9bc d9bd d9bf d9bg d9bh d9bi d9bj d9bk d9bl d9bm d9bo d9bp d9bq d9br d9bs d9bt d9bu d9bv d9bx d9by d9bz d9ca d9cb d9cc det8 mul add }if d9d 0 ne { d9d d9j d9k d9l d9n d9o d9p d9q d9r d9s d9t d9u d9w d9x d9y d9z d9aa d9ab d9ac d9ad d9af d9ag d9ah d9ai d9aj d9ak d9al d9am d9ao d9ap d9aq d9ar d9as d9at d9au d9av d9ax d9ay d9az d9ba d9bb d9bc d9bd d9be d9bg d9bh d9bi d9bj d9bk d9bl d9bm d9bn d9bp d9bq d9br d9bs d9bt d9bu d9bv d9bw d9by d9bz d9ca d9cb d9cc det8 mul sub} if d9e 0 ne { d9e d9j d9k d9l d9m d9o d9p d9q d9r d9s d9t d9u d9v d9x d9y d9z d9aa d9ab d9ac d9ad d9ae d9ag d9ah d9ai d9aj d9ak d9al d9am d9an d9ap d9aq d9ar d9as d9at d9au d9av d9aw d9ay d9az d9ba d9bb d9bc d9bd d9be d9bf d9bh d9bi d9bj d9bk d9bl d9bm d9bn d9bo d9bq d9br d9bs d9bt d9bu d9bv d9bw d9bx d9bz d9ca d9cb d9cc det8 mul add } if d9f 0 ne { d9f d9j d9k d9l d9m d9n d9p d9q d9r d9s d9t d9u d9v d9w d9y d9z d9aa d9ab d9ac d9ad d9ae d9af d9ah d9ai d9aj d9ak d9al d9am d9an d9ao d9aq d9ar d9as d9at d9au d9av d9aw d9ax d9az d9ba d9bb d9bc d9bd d9be d9bf d9bg d9bi d9bj d9bk d9bl d9bm d9bn d9bo d9bp d9br d9bs d9bt d9bu d9bv d9bw d9bx d9by d9ca d9cb d9cc det8 mul sub } if d9g 0 ne { d9g d9j d9k d9l d9m d9n d9o d9q d9r d9s d9t d9u d9v d9w d9x d9z d9aa d9ab d9ac d9ad d9ae d9af d9ag d9ai d9aj d9ak d9al d9am d9an d9ao d9ap d9ar d9as d9at d9au d9av d9aw d9ax d9ay d9ba d9bb d9bc d9bd d9be d9bf d9bg d9bh d9bj d9bk d9bl d9bm d9bn d9bo d9bp d9bq d9bs d9bt d9bu d9bv d9bw d9bx d9by d9bz d9cb d9cc det8 mul add} if d9h 0 ne { d9h d9j d9k d9l d9m d9n d9o d9p d9r d9s d9t d9u d9v d9w d9x d9y d9aa d9ab d9ac d9ad d9ae d9af d9ag d9ah d9aj d9ak d9al d9am d9an d9ao d9ap d9aq d9as d9at d9au d9av d9aw d9ax d9ay d9az d9bb d9bc d9bd d9be d9bf d9bg d9bh d9bi d9bk d9bl d9bm d9bn d9bo d9bp d9bq d9br d9bt d9bu d9bv d9bw d9bx d9by d9bz d9ca d9cc det8 mul sub} if d9i 0 ne { d9i d9j d9k d9l d9m d9n d9o d9p d9q d9s d9t d9u d9v d9w d9x d9y d9z d9ab d9ac d9ad d9ae d9af d9ag d9ah d9ai d9ak d9al d9am d9an d9ao d9ap d9aq d9ar d9at d9au d9av d9aw d9ax d9ay d9az d9ba d9bc d9bd d9be d9bf d9bg d9bh d9bi d9bj d9bl d9bm d9bn d9bo d9bp d9bq d9br d9bs d9bu d9bv d9bw d9bx d9by d9bz d9ca d9cb det8 mul add } if det9snap restore } bind def %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % COMPUTE A 10X10 FAST DETERMINANT % ================================ % /det10 accepts a b .... y z aa ab ... cb cc on the stack and % computes the resultant 10 x 10 determinant. Nine zeros are forced /det10 {save /det10snap exch def /d10cv exch def /d10cu exch def /d10ct exch def /d10cs exch def /d10cr exch def /d10cq exch def /d10cp exch def /d10co exch def /d10cn exch def /d10cm exch def /d10cl exch def /d10ck exch def /d10cj exch def /d10ci exch def /d10ch exch def /d10cg exch def /d10cf exch def /d10ce exch def /d10cd exch def /d10cc exch def /d10cb exch def /d10ca exch def /d10bz exch def /d10by exch def /d10bx exch def /d10bw exch def /d10bv exch def /d10bu exch def /d10bt exch def /d10bs exch def /d10br exch def /d10bq exch def /d10bp exch def /d10bo exch def /d10bn exch def /d10bm exch def /d10bl exch def /d10bk exch def /d10bj exch def /d10bi exch def /d10bh exch def /d10bg exch def /d10bf exch def /d10be exch def /d10bd exch def /d10bc exch def /d10bb exch def /d10ba exch def /d10az exch def /d10ay exch def /d10ax exch def /d10aw exch def /d10av exch def /d10au exch def /d10at exch def /d10as exch def /d10ar exch def /d10aq exch def /d10ap exch def /d10ao exch def /d10an exch def /d10am exch def /d10al exch def /d10ak exch def /d10aj exch def /d10ai exch def /d10ah exch def /d10ag exch def /d10af exch def /d10ae exch def /d10ad exch def /d10ac exch def /d10ab exch def /d10aa exch def /d10z exch def /d10y exch def /d10x exch def /d10w exch def /d10v exch def /d10u exch def /d10t exch def /d10s exch def /d10r exch def /d10q exch def /d10p exch def /d10o exch def /d10n exch def /d10m exch def /d10l exch def /d10k exch def /d10j exch def /d10i exch def /d10h exch def /d10g exch def /d10f exch def /d10e exch def /d10d exch def /d10c exch def /d10b exch def /d10a exch def % a b c d e f g h i j % k l m n o p q r s t % u v w x y z A B C D % E F G H I J K L M N % O P Q R S T U V W X % Y Z a b c d e f g h % i j k l m n o p q r % s t u v w x y z A B % C D E F G H I J K L % M N O P Q R S T U V d10a 0 ne { /mj d10j d10a div neg def /d10j 0 store /d10t d10k mj mul d10t add store /d10ad d10u mj mul d10ad add store /d10an d10ae mj mul d10an add store /d10ax d10ao mj mul d10ax add store /d10bh d10ay mj mul d10bh add store /d10br d10bi mj mul d10br add store /d10cb d10bs mj mul d10cb add store /d10cl d10cc mj mul d10cl add store /d10cv d10cm mj mul d10cv add store } if d10a 0 ne { /mi d10i d10a div neg def /d10i 0 store /d10s d10k mi mul d10s add store /d10ac d10u mi mul d10ac add store /d10am d10ae mi mul d10am add store /d10aw d10ao mi mul d10aw add store /d10bg d10ay mi mul d10bg add store /d10bq d10bi mi mul d10bq add store /d10ca d10bs mi mul d10ca add store /d10ck d10cc mi mul d10ck add store /d10cu d10cm mi mul d10cu add store } if d10a 0 ne { /mh d10h d10a div neg def /d10h 0 store /d10r d10k mh mul d10r add store /d10ab d10u mh mul d10ab add store /d10al d10ae mh mul d10al add store /d10av d10ao mh mul d10av add store /d10bf d10ay mh mul d10bf add store /d10bp d10bi mh mul d10bp add store /d10bz d10bs mh mul d10bz add store /d10cj d10cc mh mul d10cj add store /d10ct d10cm mh mul d10ct add store } if d10a 0 ne { /mg d10g d10a div neg def /d10g 0 store /d10q d10k mg mul d10q add store /d10aa d10u mg mul d10aa add store /d10ak d10ae mg mul d10ak add store /d10au d10ao mg mul d10au add store /d10be d10ay mg mul d10be add store /d10bo d10bi mg mul d10bo add store /d10by d10bs mg mul d10by add store /d10ci d10cc mg mul d10ci add store /d10cs d10cm mg mul d10cs add store } if d10a 0 ne { /mf d10f d10a div neg def /d10f 0 store /d10p d10k mf mul d10p add store /d10z d10u mf mul d10z add store /d10aj d10ae mf mul d10aj add store /d10at d10ao mf mul d10at add store /d10bd d10ay mf mul d10bd add store /d10bn d10bi mf mul d10bn add store /d10bx d10bs mf mul d10bx add store /d10ch d10cc mf mul d10ch add store /d10cr d10cm mf mul d10cr add store } if d10a 0 ne { /me d10e d10a div neg def /d10e 0 store /d10o d10k me mul d10o add store /d10y d10u me mul d10y add store /d10ai d10ae me mul d10ai add store /d10as d10ao me mul d10as add store /d10bc d10ay me mul d10bc add store /d10bm d10bi me mul d10bm add store /d10bw d10bs me mul d10bw add store /d10cg d10cc me mul d10cg add store /d10cq d10cm me mul d10cq add store } if d10a 0 ne { /md d10d d10a div neg def /d10d 0 store /d10n d10k md mul d10n add store /d10x d10u md mul d10x add store /d10ah d10ae md mul d10ah add store /d10ar d10ao md mul d10ar add store /d10bb d10ay md mul d10bb add store /d10bl d10bi md mul d10bl add store /d10bv d10bs md mul d10bv add store /d10cf d10cc md mul d10cf add store /d10cp d10cm md mul d10cp add store } if d10a 0 ne { /mc d10c d10a div neg def /d10c 0 store /d10m d10k mc mul d10m add store /d10w d10u mc mul d10w add store /d10ag d10ae mc mul d10ag add store /d10aq d10ao mc mul d10aq add store /d10ba d10ay mc mul d10ba add store /d10bk d10bi mc mul d10bk add store /d10bu d10bs mc mul d10bu add store /d10ce d10cc mc mul d10ce add store /d10co d10cm mc mul d10co add store } if d10a 0 ne { /mb d10b d10a div neg def /d10b 0 store /d10l d10k mb mul d10l add store /d10v d10u mb mul d10v add store /d10af d10ae mb mul d10af add store /d10ap d10ao mb mul d10ap add store /d10az d10ay mb mul d10az add store /d10bj d10bi mb mul d10bj add store /d10bt d10bs mb mul d10bt add store /d10cd d10cc mb mul d10cd add store /d10cn d10cm mb mul d10cn add store } if % a b c d e f g h i j % k l m n o p q r s t % u v w x y z A B C D % E F G H I J K L M N % O P Q R S T U V W X % Y Z a b c d e f g h % i j k l m n o p q r % s t u v w x y z A B % C D E F G H I J K L % M N O P Q R S T U V d10a 0 ne { d10a d10l d10m d10n d10o d10p d10q d10r d10s d10t d10v d10w d10x d10y d10z d10aa d10ab d10ac d10ad d10af d10ag d10ah d10ai d10aj d10ak d10al d10am d10an d10ap d10aq d10ar d10as d10at d10au d10av d10aw d10ax d10az d10ba d10bb d10bc d10bd d10be d10bf d10bg d10bh d10bj d10bk d10bl d10bm d10bn d10bo d10bp d10bq d10br d10bt d10bu d10bv d10bw d10bx d10by d10bz d10ca d10cb d10cd d10ce d10cf d10cg d10ch d10ci d10cj d10ck d10cl d10cn d10co d10cp d10cq d10cr d10cs d10ct d10cu d10cv det9 mul }{0} ifelse d10b 0 ne { d10b d10k d10m d10n d10o d10p d10q d10r d10s d10t d10u d10w d10x d10y d10z d10aa d10ab d10ac d10ad d10ae d10ag d10ah d10ai d10aj d10ak d10al d10am d10an d10ao d10aq d10ar d10as d10at d10au d10av d10aw d10ax d10ay d10ba d10bb d10bc d10bd d10be d10bf d10bg d10bh d10bi d10bk d10bl d10bm d10bn d10bo d10bp d10bq d10br d10bs d10bu d10bv d10bw d10bx d10by d10bz d10ca d10cb d10cc d10ce d10cf d10cg d10ch d10ci d10cj d10ck d10cl d10cm d10co d10cp d10cq d10cr d10cs d10ct d10cu d10cv det9 mul sub } if d10c 0 ne { d10c d10k d10l d10n d10o d10p d10q d10r d10s d10t d10u d10v d10x d10y d10z d10aa d10ab d10ac d10ad d10ae d10af d10ah d10ai d10aj d10ak d10al d10am d10an d10ao d10ap d10ar d10as d10at d10au d10av d10aw d10ax d10ay d10az d10bb d10bc d10bd d10be d10bf d10bg d10bh d10bi d10bj d10bl d10bm d10bn d10bo d10bp d10bq d10br d10bs d10bt d10bv d10bw d10bx d10by d10bz d10ca d10cb d10cc d10cd d10cf d10cg d10ch d10ci d10cj d10ck d10cl d10cm d10cn d10cp d10cq d10cr d10cs d10ct d10cu d10cv det9 mul add }if d10d 0 ne { d10d d10k d10l d10m d10o d10p d10q d10r d10s d10t d10u d10v d10w d10y d10z d10aa d10ab d10ac d10ad d10ae d10af d10ag d10ai d10aj d10ak d10al d10am d10an d10ao d10ap d10aq d10as d10at d10au d10av d10aw d10ax d10ay d10az d10ba d10bc d10bd d10be d10bf d10bg d10bh d10bi d10bj d10bk d10bm d10bn d10bo d10bp d10bq d10br d10bs d10bt d10bu d10bw d10bx d10by d10bz d10ca d10cb d10cc d10cd d10ce d10cg d10ch d10ci d10cj d10ck d10cl d10cm d10cn d10co d10cq d10cr d10cs d10ct d10cu d10cv det9 mul sub} if d10e 0 ne { d10e d10k d10l d10m d10n d10p d10q d10r d10s d10t d10u d10v d10w d10x d10z d10aa d10ab d10ac d10ad d10ae d10af d10ag d10ah d10aj d10ak d10al d10am d10an d10ao d10ap d10aq d10ar d10at d10au d10av d10aw d10ax d10ay d10az d10ba d10bb d10bd d10be d10bf d10bg d10bh d10bi d10bj d10bk d10bl d10bn d10bo d10bp d10bq d10br d10bs d10bt d10bu d10bv d10bx d10by d10bz d10ca d10cb d10cc d10cd d10ce d10cf d10ch d10ci d10cj d10ck d10cl d10cm d10cn d10co d10cp d10cr d10cs d10ct d10cu d10cv det9 mul add } if d10f 0 ne { d10f d10k d10l d10m d10n d10o d10q d10r d10s d10t d10u d10v d10w d10x d10y d10aa d10ab d10ac d10ad d10ae d10af d10ag d10ah d10ai d10ak d10al d10am d10an d10ao d10ap d10aq d10ar d10as d10au d10av d10aw d10ax d10ay d10az d10ba d10bb d10bc d10be d10bf d10bg d10bh d10bi d10bj d10bk d10bl d10bm d10bo d10bp d10bq d10br d10bs d10bt d10bu d10bv d10bw d10by d10bz d10ca d10cb d10cc d10cd d10ce d10cf d10cg d10ci d10cj d10ck d10cl d10cm d10cn d10co d10cp d10cq d10cs d10ct d10cu d10cv det9 mul sub } if d10g 0 ne { d10g d10k d10l d10m d10n d10o d10p d10r d10s d10t d10u d10v d10w d10x d10y d10z d10ab d10ac d10ad d10ae d10af d10ag d10ah d10ai d10aj d10al d10am d10an d10ao d10ap d10aq d10ar d10as d10at d10av d10aw d10ax d10ay d10az d10ba d10bb d10bc d10bd d10bf d10bg d10bh d10bi d10bj d10bk d10bl d10bm d10bn d10bp d10bq d10br d10bs d10bt d10bu d10bv d10bw d10bx d10bz d10ca d10cb d10cc d10cd d10ce d10cf d10cg d10ch d10cj d10ck d10cl d10cm d10cn d10co d10cp d10cq d10cr d10ct d10cu d10cv det9 mul add} if d10h 0 ne { d10h d10k d10l d10m d10n d10o d10p d10q d10s d10t d10u d10v d10w d10x d10y d10z d10aa d10ac d10ad d10ae d10af d10ag d10ah d10ai d10aj d10ak d10am d10an d10ao d10ap d10aq d10ar d10as d10at d10au d10aw d10ax d10ay d10az d10ba d10bb d10bc d10bd d10be d10bg d10bh d10bi d10bj d10bk d10bl d10bm d10bn d10bo d10bq d10br d10bs d10bt d10bu d10bv d10bw d10bx d10by d10ca d10cb d10cc d10cd d10ce d10cf d10cg d10ch d10ci d10ck d10cl d10cm d10cn d10co d10cp d10cq d10cr d10cs d10cu d10cv det9 mul sub} if d10i 0 ne { d10i d10k d10l d10m d10n d10o d10p d10q d10r d10t d10u d10v d10w d10x d10y d10z d10aa d10ab d10ad d10ae d10af d10ag d10ah d10ai d10aj d10ak d10al d10an d10ao d10ap d10aq d10ar d10as d10at d10au d10av d10ax d10ay d10az d10ba d10bb d10bc d10bd d10be d10bf d10bh d10bi d10bj d10bk d10bl d10bm d10bn d10bo d10bp d10br d10bs d10bt d10bu d10bv d10bw d10bx d10by d10bz d10cb d10cc d10cd d10ce d10cf d10cg d10ch d10ci d10cj d10cl d10cm d10cn d10co d10cp d10cq d10cr d10cs d10ct d10cv det9 mul add } if d10j 0 ne { d10j d10k d10l d10m d10n d10o d10p d10q d10r d10s d10u d10v d10w d10x d10y d10z d10aa d10ab d10ac d10ae d10af d10ag d10ah d10ai d10aj d10ak d10al d10am d10ao d10ap d10aq d10ar d10as d10at d10au d10av d10aw d10ay d10az d10ba d10bb d10bc d10bd d10be d10bf d10bg d10bi d10bj d10bk d10bl d10bm d10bn d10bo d10bp d10bq d10bs d10bt d10bu d10bv d10bw d10bx d10by d10bz d10ca d10cc d10cd d10ce d10cf d10cg d10ch d10ci d10cj d10ck d10cm d10cn d10co d10cp d10cq d10cr d10cs d10ct d10cu det9 mul sub } if det10snap restore } bind def %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % SOLVE A 4X4 LINEAR EQUATION % =========================== % /lineq4x4 solves an order four linear equation, given % [A B C D E F G H I J K L M N O P] and [Q R S T] % representing % % Aw + Bx + Cy + Dz = Q % Ew + Fx + Gy + Hz = R % Iw + Jx + Ky + Lz = S % Mw + Nx + Oy + Pz = T % % ...and leaving w, x, y, and z on the stack. /lineq4x4 {save /snap4x4 exch def aload pop /T4 exch def /S4 exch def /R4 exch def /Q4 exch def aload pop /P4 exch def /O4 exch def /N4 exch def /M4 exch def /L4 exch def /K4 exch def /J4 exch def /I4 exch def /H4 exch def /G4 exch def /F4 exch def /E4 exch def /D4 exch def /C4 exch def /B4 exch def /A4 exch def A4 B4 C4 D4 E4 F4 G4 H4 I4 J4 K4 L4 M4 N4 O4 P4 det4 dup Q4 B4 C4 D4 R4 F4 G4 H4 S4 J4 K4 L4 T4 N4 O4 P4 det4 exch div exch dup A4 Q4 C4 D4 E4 R4 G4 H4 I4 S4 K4 L4 M4 T4 O4 P4 det4 exch div exch dup A4 B4 Q4 D4 E4 F4 R4 H4 I4 J4 S4 L4 M4 N4 T4 P4 det4 exch div exch A4 B4 C4 Q4 E4 F4 G4 R4 I4 J4 K4 S4 M4 N4 O4 T4 det4 exch div snap4x4 restore} def % An example % 6w + 5x - 3y + 2z = 15 % -3w + x + 2y - 6z = -19 % w - 5x + 6y + 0z = 9 % 7w - 2x + 3y - 5z = -8 % % For instance % [6 5 -3 2 -3 1 2 -6 1 -5 6 0 7 -2 3 -5][15 -19 9 -8] lineq4x4 % yields w = 1, x = 2, y = 3, and z = 4. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % SOLVE A 5X5 LINEAR EQUATION % =========================== /lineq5x5 {save /snap5x5 exch def aload pop /e5 exch store /d5 exch store /c5 exch store /b5 exch store /a5 exch store aload pop /Y5 exch store /X5 exch store /W5 exch store /V5 exch store /U5 exch store /T5 exch store /S5 exch store /R5 exch store /Q5 exch store /P5 exch store /O5 exch store /N5 exch store /M5 exch store /L5 exch store /K5 exch store /J5 exch store /I5 exch store /H5 exch store /G5 exch store /F5 exch store /E5 exch store /D5 exch store /C5 exch store /B5 exch store /A5 exch store A5 B5 C5 D5 E5 F5 G5 H5 I5 J5 K5 L5 M5 N5 O5 P5 Q5 R5 S5 T5 U5 V5 W5 X5 Y5 det5 dup a5 B5 C5 D5 E5 b5 G5 H5 I5 J5 c5 L5 M5 N5 O5 d5 Q5 R5 S5 T5 e5 V5 W5 X5 Y5 det5 exch div exch dup A5 a5 C5 D5 E5 F5 b5 H5 I5 J5 K5 c5 M5 N5 O5 P5 d5 R5 S5 T5 U5 e5 W5 X5 Y5 det5 exch div exch dup A5 B5 a5 D5 E5 F5 G5 b5 I5 J5 K5 L5 c5 N5 O5 P5 Q5 d5 S5 T5 U5 V5 e5 X5 Y5 det5 exch div exch dup A5 B5 C5 a5 E5 F5 G5 H5 b5 J5 K5 L5 M5 c5 O5 P5 Q5 R5 d5 T5 U5 V5 W5 e5 Y5 det5 exch div exch A5 B5 C5 D5 a5 F5 G5 H5 I5 b5 K5 L5 M5 N5 c5 P5 Q5 R5 S5 d5 U5 V5 W5 X5 e5 det5 exch div snap5x5 restore} def % A new example % 2v + 6w + 5x - 3y + 2z = 27 % 1v - 3w + x + 2y - 6z = -24 % 3v + w - 5x + 6y + 0z = 14 % 4v + 7w - 2x + 3y - 5z = -1 % 2v + 3w - 9x + 3y - 2z = -17 % For instance % [2 6 5 -3 2 1 -3 1 2 -6 3 1 -5 6 0 4 7 -2 3 -5 % 2 3 -9 3 -2][ 27 -24 14 -1 -17] lineq5x5 % yields v = 1, w = 2, x = 3, y = 4 and z = 5. %%%%%%%%%%%% % SOLVE A 6X6 LINEAR EQUATION % =========================== /lineq6x6 {save /snap6x6 exch def aload pop /f6 exch def /e6 exch def /d6 exch def /c6 exch def /b6 exch def /a6 exch def aload pop /JJ6 exch def /II6 exch def /HH6 exch def /GG6 exch def /FF6 exch def /EE6 exch def /DD6 exch def /CC6 exch def /BB6 exch def /AA6 exch def /Z6 exch def /Y6 exch def /X6 exch def /W6 exch def /V6 exch def /U6 exch def /T6 exch def /S6 exch def /R6 exch def /Q6 exch def /P6 exch def /O6 exch def /N6 exch def /M6 exch def /L6 exch def /K6 exch def /J6 exch def /I6 exch def /H6 exch def /G6 exch def /F6 exch def /E6 exch def /D6 exch def /C6 exch def /B6 exch def /A6 exch def A6 B6 C6 D6 E6 F6 G6 H6 I6 J6 K6 L6 M6 N6 O6 P6 Q6 R6 S6 T6 U6 V6 W6 X6 Y6 Z6 AA6 BB6 CC6 DD6 EE6 FF6 GG6 HH6 II6 JJ6 det6 dup % dup report a6 B6 C6 D6 E6 F6 b6 H6 I6 J6 K6 L6 c6 N6 O6 P6 Q6 R6 d6 T6 U6 V6 W6 X6 e6 Z6 AA6 BB6 CC6 DD6 f6 FF6 GG6 HH6 II6 JJ6 det6 exch div exch dup A6 a6 C6 D6 E6 F6 G6 b6 I6 J6 K6 L6 M6 c6 O6 P6 Q6 R6 S6 d6 U6 V6 W6 X6 Y6 e6 AA6 BB6 CC6 DD6 EE6 f6 GG6 HH6 II6 JJ6 det6 exch div exch dup A6 B6 a6 D6 E6 F6 G6 H6 b6 J6 K6 L6 M6 N6 c6 P6 Q6 R6 S6 T6 d6 V6 W6 X6 Y6 Z6 e6 BB6 CC6 DD6 EE6 FF6 f6 HH6 II6 JJ6 det6 exch div exch dup A6 B6 C6 a6 E6 F6 G6 H6 I6 b6 K6 L6 M6 N6 O6 c6 Q6 R6 S6 T6 U6 d6 W6 X6 Y6 Z6 AA6 e6 CC6 DD6 EE6 FF6 GG6 f6 II6 JJ6 det6 exch div exch dup A6 B6 C6 D6 a6 F6 G6 H6 I6 J6 b6 L6 M6 N6 O6 P6 c6 R6 S6 T6 U6 V6 d6 X6 Y6 Z6 AA6 BB6 e6 DD6 EE6 FF6 GG6 HH6 f6 JJ6 det6 exch div exch A6 B6 C6 D6 E6 a6 G6 H6 I6 J6 K6 b6 M6 N6 O6 P6 Q6 c6 S6 T6 U6 V6 W6 d6 Y6 Z6 AA6 BB6 CC6 e6 EE6 FF6 GG6 HH6 II6 f6 det6 exch div snap6x6 restore} def %%%%%%%%%%%% % SOLVE A 7X7 LINEAR EQUATION % =========================== /lineq7x7 {save /snap7x7 exch def aload pop /g7 exch def /f7 exch def /e7 exch def /d7 exch def /c7 exch def /b7 exch def /a7 exch def aload pop /WW7 exch def /VV7 exch def /UU7 exch def /TT7 exch def /SS7 exch def /RR7 exch def /QQ7 exch def /PP7 exch def /OO7 exch def /NN7 exch def /MM7 exch def /LL7 exch def /KK7 exch def /JJ7 exch def /II7 exch def /HH7 exch def /GG7 exch def /FF7 exch def /EE7 exch def /DD7 exch def /CC7 exch def /BB7 exch def /AA7 exch def /Z7 exch def /Y7 exch def /X7 exch def /W7 exch def /V7 exch def /U7 exch def /T7 exch def /S7 exch def /R7 exch def /Q7 exch def /P7 exch def /O7 exch def /N7 exch def /M7 exch def /L7 exch def /K7 exch def /J7 exch def /I7 exch def /H7 exch def /G7 exch def /F7 exch def /E7 exch def /D7 exch def /C7 exch def /B7 exch def /A7 exch def A7 B7 C7 D7 E7 F7 G7 H7 I7 J7 K7 L7 M7 N7 O7 P7 Q7 R7 S7 T7 U7 V7 W7 X7 Y7 Z7 AA7 BB7 CC7 DD7 EE7 FF7 GG7 HH7 II7 JJ7 KK7 LL7 MM7 NN7 OO7 PP7 QQ7 RR7 SS7 TT7 UU7 VV7 WW7 det7 % dup a7 B7 C7 D7 E7 F7 G7 b7 I7 J7 K7 L7 M7 N7 c7 P7 Q7 R7 S7 T7 U7 d7 W7 X7 Y7 Z7 AA7 BB7 e7 DD7 EE7 FF7 GG7 HH7 II7 f7 KK7 LL7 MM7 NN7 OO7 PP7 g7 RR7 SS7 TT7 UU7 VV7 WW7 det7 exch div % exch dup dup /gotk7 exch store [ I7 J7 K7 L7 M7 N7 P7 Q7 R7 S7 T7 U7 W7 X7 Y7 Z7 AA7 BB7 DD7 EE7 FF7 GG7 HH7 II7 KK7 LL7 MM7 NN7 OO7 PP7 RR7 SS7 TT7 UU7 VV7 WW7 ] [ b7 H7 gotk7 mul sub c7 O7 gotk7 mul sub d7 V7 gotk7 mul sub e7 CC7 gotk7 mul sub f7 JJ7 gotk7 mul sub g7 QQ7 gotk7 mul sub ] lineq6x6 snap7x7 restore} def %%%%%%%%%%%%%%%% % SOLVE AN 8X8 LINEAR EQUATION % ============================ /lineq8x8 {save /snap8x8 exch def aload pop /h8 exch def /g8 exch def /f8 exch def /e8 exch def /d8 exch def /c8 exch def /b8 exch def /a8 exch def aload pop /BL8 exch def /BK8 exch def /BJ8 exch def /BI8 exch def /BH8 exch def /BG8 exch def /BF8 exch def /BE8 exch def /BD8 exch def /BC8 exch def /BB8 exch def /BA8 exch def /AZ8 exch def /AY8 exch def /AX8 exch def /AW8 exch def /AV8 exch def /AU8 exch def /AT8 exch def /AS8 exch def /AR8 exch def /AQ8 exch def /AP8 exch def /AO8 exch def /AN8 exch def /AM8 exch def /AL8 exch def /AK8 exch def /AJ8 exch def /AI8 exch def /AH8 exch def /AG8 exch def /AF8 exch def /AE8 exch def /AD8 exch def /AC8 exch def /AB8 exch def /AA8 exch def /Z8 exch def /Y8 exch def /X8 exch def /W8 exch def /V8 exch def /U8 exch def /T8 exch def /S8 exch def /R8 exch def /Q8 exch def /P8 exch def /O8 exch def /N8 exch def /M8 exch def /L8 exch def /K8 exch def /J8 exch def /I8 exch def /H8 exch def /G8 exch def /F8 exch def /E8 exch def /D8 exch def /C8 exch def /B8 exch def /A8 exch def A8 B8 C8 D8 E8 F8 G8 H8 I8 J8 K8 L8 M8 N8 O8 P8 Q8 R8 S8 T8 U8 V8 W8 X8 Y8 Z8 AA8 AB8 AC8 AD8 AE8 AF8 AG8 AH8 AI8 AJ8 AK8 AL8 AM8 AN8 AO8 AP8 AQ8 AR8 AS8 AT8 AU8 AV8 AW8 AX8 AY8 AZ8 BA8 BB8 BC8 BD8 BE8 BF8 BG8 BH8 BI8 BJ8 BK8 BL8 det8 a8 B8 C8 D8 E8 F8 G8 H8 b8 J8 K8 L8 M8 N8 O8 P8 c8 R8 S8 T8 U8 V8 W8 X8 d8 Z8 AA8 AB8 AC8 AD8 AE8 AF8 e8 AH8 AI8 AJ8 AK8 AL8 AM8 AN8 f8 AP8 AQ8 AR8 AS8 AT8 AU8 AV8 g8 AX8 AY8 AZ8 BA8 BB8 BC8 BD8 h8 BF8 BG8 BH8 BI8 BJ8 BK8 BL8 det8 exch div % exch dup dup /gotk8 exch store [ J8 K8 L8 M8 N8 O8 P8 R8 S8 T8 U8 V8 W8 X8 Z8 AA8 AB8 AC8 AD8 AE8 AF8 AH8 AI8 AJ8 AK8 AL8 AM8 AN8 AP8 AQ8 AR8 AS8 AT8 AU8 AV8 AX8 AY8 AZ8 BA8 BB8 BC8 BD8 BF8 BG8 BH8 BI8 BJ8 BK8 BL8 ] [ b8 I8 gotk8 mul sub c8 Q8 gotk8 mul sub d8 Y8 gotk8 mul sub e8 AG8 gotk8 mul sub f8 AO8 gotk8 mul sub g8 AW8 gotk8 mul sub h8 BE8 gotk8 mul sub ] lineq7x7 snap8x8 restore} def %%%%%%%%%%%%%%%% % SOLVE A 9X9 LINEAR EQUATION % =========================== /lineq9x9 {save /snap9x9 exch def aload pop /i9 exch def /h9 exch def /g9 exch def /f9 exch def /e9 exch def /d9 exch def /c9 exch def /b9 exch def /a9 exch def aload pop /CC9 exch def /CB9 exch def /CA9 exch def /BZ9 exch def /BY9 exch def /BX9 exch def /BW9 exch def /BV9 exch def /BU9 exch def /BT9 exch def /BS9 exch def /BR9 exch def /BQ9 exch def /BP9 exch def /BO9 exch def /BN9 exch def /BM9 exch def /BL9 exch def /BK9 exch def /BJ9 exch def /BI9 exch def /BH9 exch def /BG9 exch def /BF9 exch def /BE9 exch def /BD9 exch def /BC9 exch def /BB9 exch def /BA9 exch def /AZ9 exch def /AY9 exch def /AX9 exch def /AW9 exch def /AV9 exch def /AU9 exch def /AT9 exch def /AS9 exch def /AR9 exch def /AQ9 exch def /AP9 exch def /AO9 exch def /AN9 exch def /AM9 exch def /AL9 exch def /AK9 exch def /AJ9 exch def /AI9 exch def /AH9 exch def /AG9 exch def /AF9 exch def /AE9 exch def /AD9 exch def /AC9 exch def /AB9 exch def /AA9 exch def /Z9 exch def /Y9 exch def /X9 exch def /W9 exch def /V9 exch def /U9 exch def /T9 exch def /S9 exch def /R9 exch def /Q9 exch def /P9 exch def /O9 exch def /N9 exch def /M9 exch def /L9 exch def /K9 exch def /J9 exch def /I9 exch def /H9 exch def /G9 exch def /F9 exch def /E9 exch def /D9 exch def /C9 exch def /B9 exch def /A9 exch def A9 B9 C9 D9 E9 F9 G9 H9 I9 J9 K9 L9 M9 N9 O9 P9 Q9 R9 S9 T9 U9 V9 W9 X9 Y9 Z9 AA9 AB9 AC9 AD9 AE9 AF9 AG9 AH9 AI9 AJ9 AK9 AL9 AM9 AN9 AO9 AP9 AQ9 AR9 AS9 AT9 AU9 AV9 AW9 AX9 AY9 AZ9 BA9 BB9 BC9 BD9 BE9 BF9 BG9 BH9 BI9 BJ9 BK9 BL9 BM9 BN9 BO9 BP9 BQ9 BR9 BS9 BT9 BU9 BV9 BW9 BX9 BY9 BZ9 CA9 CB9 CC9 det9 a9 B9 C9 D9 E9 F9 G9 H9 I9 b9 K9 L9 M9 N9 O9 P9 Q9 R9 c9 T9 U9 V9 W9 X9 Y9 Z9 AA9 d9 AC9 AD9 AE9 AF9 AG9 AH9 AI9 AJ9 e9 AL9 AM9 AN9 AO9 AP9 AQ9 AR9 AS9 f9 AU9 AV9 AW9 AX9 AY9 AZ9 BA9 BB9 g9 BD9 BE9 BF9 BG9 BH9 BI9 BJ9 BK9 h9 BM9 BN9 BO9 BP9 BQ9 BR9 BS9 BT9 i9 BV9 BW9 BX9 BY9 BZ9 CA9 CB9 CC9 det9 exch div dup /gotk9 exch store [ K9 L9 M9 N9 O9 P9 Q9 R9 T9 U9 V9 W9 X9 Y9 Z9 AA9 AC9 AD9 AE9 AF9 AG9 AH9 AI9 AJ9 AL9 AM9 AN9 AO9 AP9 AQ9 AR9 AS9 AU9 AV9 AW9 AX9 AY9 AZ9 BA9 BB9 BD9 BE9 BF9 BG9 BH9 BI9 BJ9 BK9 BM9 BN9 BO9 BP9 BQ9 BR9 BS9 BT9 BV9 BW9 BX9 BY9 BZ9 CA9 CB9 CC9 ] [ b9 J9 gotk9 mul sub c9 S9 gotk9 mul sub d9 AB9 gotk9 mul sub e9 AK9 gotk9 mul sub f9 AT9 gotk9 mul sub g9 BC9 gotk9 mul sub h9 BL9 gotk9 mul sub i9 BU9 gotk9 mul sub ] lineq8x8 snap9x9 restore} def %%%%%%%%%%%%%%%% % SOLVE A 10X10 LINEAR EQUATION % =========================== /lineq10x10 {save /snap10x10 exch def aload pop /j10 exch def /i10 exch def /h10 exch def /g10 exch def /f10 exch def /e10 exch def /d10 exch def /c10 exch def /b10 exch def /a10 exch def aload pop /CV10 exch def /CU10 exch def /CT10 exch def /CS10 exch def /CR10 exch def /CQ10 exch def /CP10 exch def /CO10 exch def /CN10 exch def /CM10 exch def /CL10 exch def /CK10 exch def /CJ10 exch def /CI10 exch def /CH10 exch def /CG10 exch def /CF10 exch def /CE10 exch def /CD10 exch def /CC10 exch def /CB10 exch def /CA10 exch def /BZ10 exch def /BY10 exch def /BX10 exch def /BW10 exch def /BV10 exch def /BU10 exch def /BT10 exch def /BS10 exch def /BR10 exch def /BQ10 exch def /BP10 exch def /BO10 exch def /BN10 exch def /BM10 exch def /BL10 exch def /BK10 exch def /BJ10 exch def /BI10 exch def /BH10 exch def /BG10 exch def /BF10 exch def /BE10 exch def /BD10 exch def /BC10 exch def /BB10 exch def /BA10 exch def /AZ10 exch def /AY10 exch def /AX10 exch def /AW10 exch def /AV10 exch def /AU10 exch def /AT10 exch def /AS10 exch def /AR10 exch def /AQ10 exch def /AP10 exch def /AO10 exch def /AN10 exch def /AM10 exch def /AL10 exch def /AK10 exch def /AJ10 exch def /AI10 exch def /AH10 exch def /AG10 exch def /AF10 exch def /AE10 exch def /AD10 exch def /AC10 exch def /AB10 exch def /AA10 exch def /Z10 exch def /Y10 exch def /X10 exch def /W10 exch def /V10 exch def /U10 exch def /T10 exch def /S10 exch def /R10 exch def /Q10 exch def /P10 exch def /O10 exch def /N10 exch def /M10 exch def /L10 exch def /K10 exch def /J10 exch def /I10 exch def /H10 exch def /G10 exch def /F10 exch def /E10 exch def /D10 exch def /C10 exch def /B10 exch def /A10 exch def A10 B10 C10 D10 E10 F10 G10 H10 I10 J10 K10 L10 M10 N10 O10 P10 Q10 R10 S10 T10 U10 V10 W10 X10 Y10 Z10 AA10 AB10 AC10 AD10 AE10 AF10 AG10 AH10 AI10 AJ10 AK10 AL10 AM10 AN10 AO10 AP10 AQ10 AR10 AS10 AT10 AU10 AV10 AW10 AX10 AY10 AZ10 BA10 BB10 BC10 BD10 BE10 BF10 BG10 BH10 BI10 BJ10 BK10 BL10 BM10 BN10 BO10 BP10 BQ10 BR10 BS10 BT10 BU10 BV10 BW10 BX10 BY10 BZ10 CA10 CB10 CC10 CD10 CE10 CF10 CG10 CH10 CI10 CJ10 CK10 CL10 CM10 CN10 CO10 CP10 CQ10 CR10 CS10 CT10 CU10 CV10 det10 a10 B10 C10 D10 E10 F10 G10 H10 I10 J10 b10 L10 M10 N10 O10 P10 Q10 R10 S10 T10 c10 V10 W10 X10 Y10 Z10 AA10 AB10 AC10 AD10 d10 AF10 AG10 AH10 AI10 AJ10 AK10 AL10 AM10 AN10 e10 AP10 AQ10 AR10 AS10 AT10 AU10 AV10 AW10 AX10 f10 AZ10 BA10 BB10 BC10 BD10 BE10 BF10 BG10 BH10 g10 BJ10 BK10 BL10 BM10 BN10 BO10 BP10 BQ10 BR10 h10 BT10 BU10 BV10 BW10 BX10 BY10 BZ10 CA10 CB10 i10 CD10 CE10 CF10 CG10 CH10 CI10 CJ10 CK10 CL10 j10 CN10 CO10 CP10 CQ10 CR10 CS10 CT10 CU10 CV10 det10 exch div dup /gotk10 exch store [ L10 M10 N10 O10 P10 Q10 R10 S10 T10 V10 W10 X10 Y10 Z10 AA10 AB10 AC10 AD10 AF10 AG10 AH10 AI10 AJ10 AK10 AL10 AM10 AN10 AP10 AQ10 AR10 AS10 AT10 AU10 AV10 AW10 AX10 AZ10 BA10 BB10 BC10 BD10 BE10 BF10 BG10 BH10 BJ10 BK10 BL10 BM10 BN10 BO10 BP10 BQ10 BR10 BT10 BU10 BV10 BW10 BX10 BY10 BZ10 CA10 CB10 CD10 CE10 CF10 CG10 CH10 CI10 CJ10 CK10 CL10 CN10 CO10 CP10 CQ10 CR10 CS10 CT10 CU10 CV10 ] [ b10 K10 gotk10 mul sub c10 U10 gotk10 mul sub d10 AE10 gotk10 mul sub e10 AO10 gotk10 mul sub f10 AY10 gotk10 mul sub g10 BI10 gotk10 mul sub h10 BS10 gotk10 mul sub i10 CC10 gotk10 mul sub j10 CM10 gotk10 mul sub ] lineq9x9 snap10x10 restore} def %%%%%%%%%%%%%%%%%%%%% SOME DEMOS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Linear phase digital filters make fine demos of using determinants % to solve high variable linear equations. % First some log plotting utilities... % /makeloggraph does a 0 to 50 frequency 0 to -60 decibel background plot. % Note that the position and scale (nominal 10x) has to be pretranslated. /makealoggraph {0.012 setlinewidth 1 setlinecap 1 setlinejoin 1 1 50 { /posn exch store posn log 30 mul 0 moveto 0 60 rlineto stroke} for gsave 61 {0 0 moveto 50 log 30 mul 0 lineto stroke 0 1 translate} repeat grestore 0.10 setlinewidth gsave 7 {0 0 moveto 50 log 30 mul 0 lineto stroke 0 10 translate} repeat grestore 1 log 30 mul 0 moveto 0 60 rlineto stroke 10 log 30 mul 0 moveto 0 60 rlineto stroke 20 log 30 mul 0 moveto 0 60 rlineto stroke 30 log 30 mul 0 moveto 0 60 rlineto stroke 40 log 30 mul 0 moveto 0 60 rlineto stroke 50 log 30 mul 0 moveto 0 60 rlineto stroke} def % This plots a point on the graph /mdot {newpath 0.150 0 360 arc fill } def %% % alternate to return recordable data to host %% /mdot { 10 mul cvi 10 div == 10 mul cvi 10 div ==} def % Here is a second order Butterworth lowpass filter with a 2f cutoff freq. % This is useful for comparisons. You can also replace it with your desired % target filter response. /Butterworth2x {1 0.5 50 {/freq exch store freq 4 div dup mul dup mul 1 add sqrt 1 exch div freq log 30 mul exch log 20 mul 60 add mdot} for} def % Make a graph and a Butterworth plot.. (uncomment to view) 100 100 translate 8 dup scale makealoggraph %% Butterworth2x %%%%%%%%%%%%%%%%%% LINEAR PHASE FILTERS %%%%%%%%%%%%%%%%%%%%%% % In a linear phase digital filter, paired samples are taken in a time % interval t. Multiples of t combined with proper coefficients generate the % desired response. Here are some plotting utilities... % /plotlinephase3 plots the response of a three coefficient filter of % coeffB(priorsamp1) + coeffA + coeffB(postsamp1) /sampfreq 72 def % clock sampling frequency /plotfreq 0.1 def % plot dot spacing /sampphase {360 sampfreq div} def % unity frequency sample phase /scaleit {sampphase mul mul cos} def % convenience scaling operator /plotit {abs dup 0 eq {pop 0.0000001} % convenience plotting operator if log 20 mul 60 add freq log 30 mul exch mdot} def /plotlinphase3 {1 plotfreq 50 {/freq exch store freq 1 scaleit coeffB mul 2 mul coeffA add plotit} for} def % use example... (uncomment to view) %% /coeffA 0.500 def /coeffB 0.250 def plotlinphase3 % /plotlinephase5 plots the response of a five coefficient filter of % coeffC(priorsamp2) + coeffB(priorsamp1) + coeffA + % coeffB(postsamp1) + coeffC(postsamp2) /plotlinphase5 {1 plotfreq 50 {/freq exch store freq 2 scaleit coeffC mul 2 mul freq 1 scaleit coeffB mul 2 mul add coeffA add plotit} for} def % and the rest of them... /plotlinphase7 {1 plotfreq 50 {/freq exch store freq 3 scaleit coeffD mul 2 mul freq 2 scaleit coeffC mul 2 mul add freq 1 scaleit coeffB mul 2 mul add coeffA add plotit} for} def /plotlinphase9 {1 plotfreq 50 {/freq exch store freq 4 scaleit coeffE mul 2 mul freq 3 scaleit coeffD mul 2 mul add freq 2 scaleit coeffC mul 2 mul add freq 1 scaleit coeffB mul 2 mul add coeffA add plotit} for} def /plotlinphase11 {1 plotfreq 50 {/freq exch store freq 5 scaleit coeffF mul 2 mul freq 4 scaleit coeffE mul 2 mul add freq 3 scaleit coeffD mul 2 mul add freq 2 scaleit coeffC mul 2 mul add freq 1 scaleit coeffB mul 2 mul add coeffA add plotit} for} def /plotlinphase13 {1 plotfreq 50 {/freq exch store freq 6 scaleit coeffG mul 2 mul freq 5 scaleit coeffF mul 2 mul add freq 4 scaleit coeffE mul 2 mul add freq 3 scaleit coeffD mul 2 mul add freq 2 scaleit coeffC mul 2 mul add freq 1 scaleit coeffB mul 2 mul add coeffA add plotit} for} def /plotlinphase15 {1 plotfreq 50 {/freq exch store freq 7 scaleit coeffH mul 2 mul freq 6 scaleit coeffG mul 2 mul add freq 5 scaleit coeffF mul 2 mul add freq 4 scaleit coeffE mul 2 mul add freq 3 scaleit coeffD mul 2 mul add freq 2 scaleit coeffC mul 2 mul add freq 1 scaleit coeffB mul 2 mul add coeffA add plotit} for} def /plotlinphase17 {1 plotfreq 50 {/freq exch store freq 8 scaleit coeffI mul 2 mul freq 7 scaleit coeffH mul 2 mul add freq 6 scaleit coeffG mul 2 mul add freq 5 scaleit coeffF mul 2 mul add freq 4 scaleit coeffE mul 2 mul add freq 3 scaleit coeffD mul 2 mul add freq 2 scaleit coeffC mul 2 mul add freq 1 scaleit coeffB mul 2 mul add coeffA add plotit} for} def /plotlinphase19 {1 plotfreq 50 {/freq exch store freq 9 scaleit coeffJ mul 2 mul freq 8 scaleit coeffI mul 2 mul add freq 7 scaleit coeffH mul 2 mul add freq 6 scaleit coeffG mul 2 mul add freq 5 scaleit coeffF mul 2 mul add freq 4 scaleit coeffE mul 2 mul add freq 3 scaleit coeffD mul 2 mul add freq 2 scaleit coeffC mul 2 mul add freq 1 scaleit coeffB mul 2 mul add coeffA add plotit} for} def % FINDING COEFFICIENTS % ==================== % Linear algebraic equations let us find the needed coefficients for % any order filter. Consider a four coefficient (seven sample) filter... /freqs [1 8 20 30] def % set your frequencies here Don't use zero! /amps [1 0.707 0 0] def % set your amplitudes here /slp2 {sampphase mul mul cos 2 mul} def % convenience operator /findcoeffs4 { /f1 freqs 0 get store /f2 freqs 1 get store /f3 freqs 2 get store /f4 freqs 3 get store [ f1 3 slp2 f1 2 slp2 f1 1 slp2 1 % set up frequency equations f2 3 slp2 f2 2 slp2 f2 1 slp2 1 f3 3 slp2 f3 2 slp2 f3 1 slp2 1 f4 3 slp2 f4 2 slp2 f4 1 slp2 1 ] amps % get desired amplitudes lineq4x4 % solve linear 4x4 equation /coeffA exch store % grab the coefficients /coeffB exch store /coeffC exch store /coeffD exch store plotlinphase7 % and flaunt them } bind def % a demo... (uncomment to display) %% /freqs [1 8 20 30] def % set your frequencies here Don't use zero! %% /amps [1 0.707 0 0] def % set your amplitudes here %% findcoeffs4 % find and plot your response %%%%%%%%%%%%%%%%%% % Five coefficients, nine samples... /findcoeffs5 { /f1 freqs 0 get store /f2 freqs 1 get store /f3 freqs 2 get store /f4 freqs 3 get store /f5 freqs 4 get store [ f1 4 slp2 f1 3 slp2 f1 2 slp2 f1 1 slp2 1 % set up frequency equations f2 4 slp2 f2 3 slp2 f2 2 slp2 f2 1 slp2 1 f3 4 slp2 f3 3 slp2 f3 2 slp2 f3 1 slp2 1 f4 4 slp2 f4 3 slp2 f4 2 slp2 f4 1 slp2 1 f5 4 slp2 f5 3 slp2 f5 2 slp2 f5 1 slp2 1 ] amps % get desired amplitudes lineq5x5 % solve linear 4x4 equation /coeffA exch store % grab the coefficients /coeffB exch store /coeffC exch store /coeffD exch store /coeffE exch store plotlinphase9 % and flaunt them } bind def % demo... (uncomment to display) %% /freqs [1 8 16 22 30] def % set your frequencies here Don't use zero! %% /amps [1 0.707 0 0 0] def % set your amplitudes here %% findcoeffs5 % find and plot your response %%%%%%%%%%%%%%%%%% % Six coefficients, eleven samples... /findcoeffs6 { /f1 freqs 0 get store /f2 freqs 1 get store /f3 freqs 2 get store /f4 freqs 3 get store /f5 freqs 4 get store /f6 freqs 5 get store [ f1 5 slp2 f1 4 slp2 f1 3 slp2 f1 2 slp2 f1 1 slp2 1 % set up frequency equations f2 5 slp2 f2 4 slp2 f2 3 slp2 f2 2 slp2 f2 1 slp2 1 f3 5 slp2 f3 4 slp2 f3 3 slp2 f3 2 slp2 f3 1 slp2 1 f4 5 slp2 f4 4 slp2 f4 3 slp2 f4 2 slp2 f4 1 slp2 1 f5 5 slp2 f5 4 slp2 f5 3 slp2 f5 2 slp2 f5 1 slp2 1 f6 5 slp2 f6 4 slp2 f6 3 slp2 f6 2 slp2 f6 1 slp2 1 ] amps % get desired amplitudes lineq6x6 % solve linear 4x4 equation /coeffA exch store % grab the coefficients /coeffB exch store /coeffC exch store /coeffD exch store /coeffE exch store /coeffF exch store plotlinphase11 % and flaunt them } bind def % demo... (uncomment to display) %% /freqs [1 8 16 20 24 30] def % set your frequencies here Don't use zero! %% /amps [1 0.707 0 0 0 0] def % set your amplitudes here %% findcoeffs6 %%%%%%%%%%%%%%%%%% % Seven coefficients, thirteen samples... /findcoeffs7 { /f1 freqs 0 get store /f2 freqs 1 get store /f3 freqs 2 get store /f4 freqs 3 get store /f5 freqs 4 get store /f6 freqs 5 get store /f7 freqs 6 get store [ f1 6 slp2 f1 5 slp2 f1 4 slp2 f1 3 slp2 f1 2 slp2 f1 1 slp2 1 f2 6 slp2 f2 5 slp2 f2 4 slp2 f2 3 slp2 f2 2 slp2 f2 1 slp2 1 f3 6 slp2 f3 5 slp2 f3 4 slp2 f3 3 slp2 f3 2 slp2 f3 1 slp2 1 f4 6 slp2 f4 5 slp2 f4 4 slp2 f4 3 slp2 f4 2 slp2 f4 1 slp2 1 f5 6 slp2 f5 5 slp2 f5 4 slp2 f5 3 slp2 f5 2 slp2 f5 1 slp2 1 f6 6 slp2 f6 5 slp2 f6 4 slp2 f6 3 slp2 f6 2 slp2 f6 1 slp2 1 f7 6 slp2 f7 5 slp2 f7 4 slp2 f7 3 slp2 f7 2 slp2 f7 1 slp2 1 ] amps % get desired amplitudes lineq7x7 % solve linear 4x4 equation /coeffA exch store % grab the coefficients /coeffB exch store /coeffC exch store /coeffD exch store /coeffE exch store /coeffF exch store /coeffG exch store plotlinphase13 % and flaunt them } def %% /freqs [1 8 16 20 24 28 32] def % set your frequencies here %% /amps [1 0.707 0 0 0 0 0] def % set your amplitudes here %% findcoeffs7 %% /freqs [1 6 14 18 22 26 30] def % set your frequencies here %% /amps [1 0.707 0.01 0.01 0.01 0.01 0.01 ] def % set your amplitudes here %% findcoeffs7 %%%%%%%%%%%%%%%%%% % Eight coefficients, fifteen samples... /findcoeffs8 { /f1 freqs 0 get store /f2 freqs 1 get store /f3 freqs 2 get store /f4 freqs 3 get store /f5 freqs 4 get store /f6 freqs 5 get store /f7 freqs 6 get store /f8 freqs 7 get store [ f1 7 slp2 f1 6 slp2 f1 5 slp2 f1 4 slp2 f1 3 slp2 f1 2 slp2 f1 1 slp2 1 f2 7 slp2 f2 6 slp2 f2 5 slp2 f2 4 slp2 f2 3 slp2 f2 2 slp2 f2 1 slp2 1 f3 7 slp2 f3 6 slp2 f3 5 slp2 f3 4 slp2 f3 3 slp2 f3 2 slp2 f3 1 slp2 1 f4 7 slp2 f4 6 slp2 f4 5 slp2 f4 4 slp2 f4 3 slp2 f4 2 slp2 f4 1 slp2 1 f5 7 slp2 f5 6 slp2 f5 5 slp2 f5 4 slp2 f5 3 slp2 f5 2 slp2 f5 1 slp2 1 f6 7 slp2 f6 6 slp2 f6 5 slp2 f6 4 slp2 f6 3 slp2 f6 2 slp2 f6 1 slp2 1 f7 7 slp2 f7 6 slp2 f7 5 slp2 f7 4 slp2 f7 3 slp2 f7 2 slp2 f7 1 slp2 1 f8 7 slp2 f8 6 slp2 f8 5 slp2 f8 4 slp2 f8 3 slp2 f8 2 slp2 f8 1 slp2 1 ] amps % get desired amplitudes lineq8x8 % solve linear 4x4 equation /coeffA exch store % grab the coefficients /coeffB exch store /coeffC exch store /coeffD exch store /coeffE exch store /coeffF exch store /coeffG exch store /coeffH exch store plotlinphase15 % and flaunt them } def %% /freqs [1 6 14 18 22 26 30 34] def % set your frequencies here %% /amps [1 0.707 0.003 0.003 0.003 0.003 0.003 0.003] def % set amps %% findcoeffs8 %%%%%%%%%%%%%%%%%% % Nine coefficients, seventeen samples... /findcoeffs9 { /f1 freqs 0 get store /f2 freqs 1 get store /f3 freqs 2 get store /f4 freqs 3 get store /f5 freqs 4 get store /f6 freqs 5 get store /f7 freqs 6 get store /f8 freqs 7 get store /f9 freqs 8 get store [ f1 8 slp2 f1 7 slp2 f1 6 slp2 f1 5 slp2 f1 4 slp2 f1 3 slp2 f1 2 slp2 f1 1 slp2 1 f2 8 slp2 f2 7 slp2 f2 6 slp2 f2 5 slp2 f2 4 slp2 f2 3 slp2 f2 2 slp2 f2 1 slp2 1 f3 8 slp2 f3 7 slp2 f3 6 slp2 f3 5 slp2 f3 4 slp2 f3 3 slp2 f3 2 slp2 f3 1 slp2 1 f4 8 slp2 f4 7 slp2 f4 6 slp2 f4 5 slp2 f4 4 slp2 f4 3 slp2 f4 2 slp2 f4 1 slp2 1 f5 8 slp2 f5 7 slp2 f5 6 slp2 f5 5 slp2 f5 4 slp2 f5 3 slp2 f5 2 slp2 f5 1 slp2 1 f6 8 slp2 f6 7 slp2 f6 6 slp2 f6 5 slp2 f6 4 slp2 f6 3 slp2 f6 2 slp2 f6 1 slp2 1 f7 8 slp2 f7 7 slp2 f7 6 slp2 f7 5 slp2 f7 4 slp2 f7 3 slp2 f7 2 slp2 f7 1 slp2 1 f8 8 slp2 f8 7 slp2 f8 6 slp2 f8 5 slp2 f8 4 slp2 f8 3 slp2 f8 2 slp8 f8 1 slp2 1 f9 8 slp2 f9 7 slp2 f9 6 slp2 f9 5 slp2 f9 4 slp2 f9 3 slp2 f9 2 slp8 f9 1 slp2 1 ] amps % get desired amplitudes lineq9x9 % solve linear 9x9 equation /coeffA exch store % grab the coefficients /coeffB exch store /coeffC exch store /coeffD exch store /coeffE exch store /coeffF exch store /coeffG exch store /coeffH exch store /coeffI exch store plotlinphase17 % and flaunt them } def %%%%%%%%%%%%%%%%%% %% /freqs [1 6 16.2 %% 19.5 22.5 25.5 28.5 31.5 34.5] def % set your frequencies here %% /amps [1 0.707 0.003 0.003 0.003 0.003 0.003 0.003 0.003 ] def % set amps %% findcoeffs9 %%%%%%%%%%%%%%%%%% % Nine coefficients, seventeen samples... /findcoeffs9 { /f1 freqs 0 get store /f2 freqs 1 get store /f3 freqs 2 get store /f4 freqs 3 get store /f5 freqs 4 get store /f6 freqs 5 get store /f7 freqs 6 get store /f8 freqs 7 get store /f9 freqs 8 get store /f10 freqs 9 get store [ f1 9 slp2 f1 8 slp2 f1 7 slp2 f1 6 slp2 f1 5 slp2 f1 4 slp2 f1 3 slp2 f1 2 slp2 f1 1 slp2 1 f2 9 slp2 f2 8 slp2 f2 7 slp2 f2 6 slp2 f2 5 slp2 f2 4 slp2 f2 3 slp2 f2 2 slp2 f2 1 slp2 1 f3 9 slp2 f3 8 slp2 f3 7 slp2 f3 6 slp2 f3 5 slp2 f3 4 slp2 f3 3 slp2 f3 2 slp2 f3 1 slp2 1 f4 9 slp2 f4 8 slp2 f4 7 slp2 f4 6 slp2 f4 5 slp2 f4 4 slp2 f4 3 slp2 f4 2 slp2 f4 1 slp2 1 f5 9 slp2 f5 8 slp2 f5 7 slp2 f5 6 slp2 f5 5 slp2 f5 4 slp2 f5 3 slp2 f5 2 slp2 f5 1 slp2 1 f6 9 slp2 f6 8 slp2 f6 7 slp2 f6 6 slp2 f6 5 slp2 f6 4 slp2 f6 3 slp2 f6 2 slp2 f6 1 slp2 1 f7 9 slp2 f7 8 slp2 f7 7 slp2 f7 6 slp2 f7 5 slp2 f7 4 slp2 f7 3 slp2 f7 2 slp2 f7 1 slp2 1 f8 9 slp2 f8 8 slp2 f8 7 slp2 f8 6 slp2 f8 5 slp2 f8 4 slp2 f8 3 slp2 f8 2 slp2 f8 1 slp2 1 f9 9 slp2 f9 8 slp2 f9 7 slp2 f9 6 slp2 f9 5 slp2 f9 4 slp2 f9 3 slp2 f9 2 slp2 f9 1 slp2 1 f10 9 slp2 f10 8 slp2 f10 7 slp2 f10 6 slp2 f10 5 slp2 f10 4 slp2 f10 3 slp2 f10 2 slp2 f10 1 slp2 1 ] amps % get desired amplitudes lineq10x10 % solve linear 10x10 equation /coeffA exch store % grab the coefficients /coeffB exch store /coeffC exch store /coeffD exch store /coeffE exch store /coeffF exch store /coeffG exch store /coeffH exch store /coeffI exch store /coeffJ exch store plotlinphase19 % and flaunt them } def %%%%%%%%%%%%%%%%%% /freqs [1 6 13.5 0.8 add 16.5 0.5 add 19.5 0.2 add 22.5 25.5 28.5 31.5 34.5] def % set your frequencies here % good set with bounce % /freqs [1 6 14 16 18.7 21.9 24.75 28.35 31.45 33.2 ] def % /amps [1 0.707 0.004 0.004 0.004 0.004 0.004 0.004 0.004 0.004 ] def % set amps % way down in the 70's! % /freqs [1 6 14 15.7 18.5 21.3 24.4 27.9 31.1 34.2] def % /amps [1 .707 0 0 0 0 0 0 0 0] def % tight but only -33 db %% /freqs [1 6 10 12.25 15.5 19.1 22.8 26.8 30.4 34.1] def %% /amps [1 .707 0 0 0 0 0 0 0 0] def % nice at -41 db down /freqs [1 6 11 12.95 15.92 19.32 22.95 26.6 30.35 34.1] def /amps [1 .707 0 0 0 0 0 0 0 0] def findcoeffs9 %%%%%%%%%%%%%%%%%% showpage quit