home *** CD-ROM | disk | FTP | other *** search
/ The C Users' Group Library 1994 August / wc-cdrom-cusersgrouplibrary-1994-08.iso / vol_200 / 209_01 / simpmain.c < prev    next >
Text File  |  1990-03-04  |  10KB  |  354 lines
  1. /* SIMPMAIN.C    VERS:- 01.00  DATE:- 09/26/86  TIME:- 09:39:29 PM */
  2. /* 
  3. Description:
  4.  
  5. main program controlling input, output & simplex fitting     = main()
  6.  
  7.  
  8. By J.A« Rupley¼ Tucson¼ Arizona
  9. Coded for ECO C compiler, version 3.40
  10. */
  11.  
  12. /*
  13. Somσ comment≤ oε thσ constructioε oµ SIMPMAIN anΣ it≤ ì
  14. supporting modules.
  15.  
  16.  
  17.      Abou⌠ 4╦ oµ memor∙ arσ reserveΣ t∩ allo≈ expansioε oµ ì
  18. <main()╛ and <func()>¼ fo≥ morσ elaboratσ output¼ additional 
  19. functions¼ etc.
  20.  
  21.  
  22.      Thσ maximuφ numbe≥ oµ parameters¼ NPARM¼ i≤ se⌠ a⌠ 10; iµ ì
  23. morσ arσ needed¼ al∞ routine≤ (SIMPMAIN¼ SIMPLIB0¼ etc⌐ mus⌠ bσ ì
  24. recompileΣ witΦ thσ ne≈ valuσ oµ NPARM.
  25.  
  26.  
  27.      T∩ makσ morσ readablσ thσ codinτ iε <func()╛ oµ thσ mode∞ ì
  28. equatioε t∩ bσ fi⌠ t∩ thσ data║ 
  29.  
  30. (1⌐ usσ mnemoniπ membe≥ name≤ iε declarinτ <struc⌠ dat╛ iε ì
  31. XXXXFITn╗ 
  32.  
  33. (2⌐ declarσ ß dumm∙ structure¼ <struc⌠ pnamestruct>¼ tha⌠ i≤ ì
  34. entirel∙ equivalen⌠ t∩ thσ structurσ tha⌠ hold≤ thσ paramete≥ ì
  35. values¼ <pstruct>¼ bu⌠ tha⌠ ha≤ mnemoniπ membe≥ names╗ thσ ì
  36. mnemoniπ dumm∙ structurσ theε caε bσ useΣ with thσ <pstruct╛ ì
  37. addres≤ passeΣ a≤ ß paramete≥ t∩ <func()>.
  38.  
  39.  
  40.      Thσ DEFINITIO╬ oµ thσ aggregatσ <data╛ anΣ thσ function≤ ì
  41. <use_mess()>¼ <file()>¼ anΣ <read_data()╛ arσ iε ß separatσ file¼ ì
  42. SIMPLIB1╗ thi≤ i≤ t∩ t∩ allo≈ expansioε oµ thσ aggregatσ <data╛ ì
  43. b∙ overwritinτ mos⌠ oµ SIMPLIB1╗ thσ SIMPLIB▒ routine≤ arσ ì
  44. entereΣ onl∙ once¼ a⌠ thσ star⌠ oµ execution.
  45.  
  46.  
  47.      thσ DECLARATIO╬ oµ <struc⌠ dat╛ accordinτ t∩ thσ ì
  48. requirement≤ oµ thσ mode∞ describeΣ b∙ <func()╛ i≤ giveε iε ì
  49. XXXXFITn╗ <func()>¼ etc« referencσ thσ aggregatσ <data╛ a≤ ì
  50. externa∞ t∩ XXXXFITn¼ bu⌠ througΦ thσ structurσ <dat>¼ declareΣ ì
  51. locall∙ iε XXXXFITn witΦ mnemoniπ membe≥ name≤ suitablσ fo≥ usσ ì
  52. iε thσ codinτ oµ <func()>¼ etc.
  53.  
  54.  
  55.      Thσ inten⌠ i≤ t∩ generalizσ thσ reaΣ oµ thσ datß filσ anΣ ì
  56. thσ allocatioε oµ datß storagσ (iε SIMPLIB1)¼ whilσ retaininτ ì
  57. flexibilit∙ iε thσ declaratioε oµ <struc⌠ da⌠ data╛ iε XXXXFITn╗ ì
  58. thσ followinτ comment≤ bea≥ upoε thi≤ arrangement:
  59.  
  60.  
  61.      Thσ loadinτ oµ value≤ int∩ thσ aggregatσ <data╛ i≤ donσ iε ì
  62. thσ SIMPLIB▒ modulσ <read_data()> b∙ usσ of:
  63.  
  64. (1⌐ ß generalizeΣ ("dummy"⌐ structurσ fo≥ <data>;
  65.  
  66. (2⌐ ß reaΣ loo≡ tha⌠ move≤ successivσ doublσ value≤ froφ thσ ì
  67. asciΘ datß filσ int∩ thσ storagσ a⌠ anΣ abovσ <data[0]>¼ withou⌠ ì
  68. referencinτ thσ element≤ oµ <data╛ b∙ structurσ membe≥ o≥ index.
  69.  
  70.  
  71.      Thσ "usefuló declaratioε oµ thσ structurσ fo≥ <data╛ i≤ iε ì
  72. XXXXFITn¼ wherσ i⌠ i≤ referenceΣ b∙ <func()╛ anΣ <fdatprint()>╗ ì
  73. <struc⌠ dat╛ mus⌠ bσ changeΣ t∩ accorΣ witΦ thσ requirement≤ oµ ì
  74. <func()╛ anΣ <fdatprint()>╗ al∞ member≤ oµ <struc⌠ dat╛ MUS╘ bσ ì
  75. oµ typσ double.
  76.  
  77.  
  78.      Changσ iε thσ mode∞ beinτ fi⌠ shoulΣ no⌠ requirσ recodinτ ì
  79. anΣ recompilatioε oµ <read_data()╛ o≥ oµ an∙ othe≥ routine≤ excep⌠ ì
  80. thosσ oµ XXXXFITn╗ oµ course¼ changσ iε thσ mode∞ require≤ changσ ì
  81. oµ <func()>¼ <fdatprint()>¼ anΣ of thσ declaratioεs oµ <struc⌠ dat╛ ì
  82. anΣ <struc⌠ pnamestruct╛ iε XXXXFITn.
  83.  
  84.  
  85.      Thσ sizσ oµ thσ <data╛ aggregatσ i≤ limiteΣ b∙ (1⌐ thσ sizσ ì
  86. oµ freσ memor∙ anΣ (2⌐ thσ sizσ oµ SIMPLIB1¼ mos⌠ oµ whicΦ caε bσ ì
  87. overwritteε b∙ datß records╗ fo≥ thi≤ versioε oµ thσ program¼ ì
  88. SIMPLIB▒ correspond≤ t∩ abou⌠ 60░ doublσ values¼ anΣ unuseΣ ì
  89. memor∙ t∩ abou⌠ 60░ doublσ values╗ overwritinτ oµ thσ codσ oµ ì
  90. SIMPLIB▒ ma∙ no⌠ bσ alloweΣ b∙ somσ compilers.
  91.  
  92.  
  93.      Fo≥ thσ six-membe≥ structurσ <dat╛ useΣ iε thi≤ versioε oµ ì
  94. thσ program¼ thσ maximuφ numbe≥ oµ datß point≤ i≤ morσ thaε 10░ ì
  95. (60░ doublσ values)¼ expandablσ t∩ morσ thaε 20░ (120░ doublσ ì
  96. values⌐ iµ SIMPLIB▒ i≤ recompileΣ witΦ aε increasσ oµ NDATA╗ ì
  97. NDAT┴ i≤ currentl∙ se⌠ a⌠ 35░ doublσ values╗ increasσ oµ NDAT┴ oµ ì
  98. coursσ decrease≤ thσ amoun⌠ oµ memor∙ availablσ fo≥ expansioε oµ ì
  99. thσ codσ oµ <main()>¼ <func()>¼ etc.
  100. */
  101.  
  102. /* page eject */
  103.  
  104. #include <stdio.h>
  105. #include <ctrlcnst.h>
  106.  
  107. #define NPARM        10    /* do NOT change this define */
  108.  
  109.  
  110.  
  111.                 /* STRUCTURES */
  112.  
  113.                 /* do NOT change any structure */
  114.  
  115. struct pstruct {
  116.     double val ;
  117.     double parm[NPARM] ;
  118. } ;
  119.  
  120. struct qstruct {
  121.     int parmndx[NPARM] ;
  122.     double q[NPARM] ;
  123.     double yplus[NPARM] ;
  124.     double yminus[NPARM] ;
  125.     double a[NPARM] ;
  126.     double bdiag[NPARM] ;
  127.     double inv_bdiag[NPARM] ;
  128.     double std_dev[NPARM] ;
  129. } ;
  130.  
  131. /* page eject */
  132.  
  133. /* MAIN
  134.     MAIN PROGRAM FOR CONTROL OF:
  135.     DATA INPUT FROM DISK FILE = ARGUMENT 1 OF COMMAND LINE
  136.     SIMPLEX FITTING
  137.     QUADRATIC FIT FOR EXTRACTION OF STANDARD DEVIATIONS
  138.     OUTPUT TO CONSOLE AND, IF SPECIFIED IN OPTIONAL ARGUMENT 2 OF
  139.         COMMAND LINE, TO PRINTER OR DISKFILE
  140. */
  141.  
  142. /* 
  143. at entry, one must initialize the following, by a call to <read_data()>:
  144. double exit_test, quad_test
  145. int prt_cycle, maxquad_skip
  146. int iter, maxiter, nparm, nvert, nfree, ndatval, ndata
  147. char title[80]
  148. struc⌠ pstruc⌠ p[nvert▌ ╜ thσ startinτ simplex
  149. struct dat data[ndata]    = the data array, used in <func()> and <fdatprint()>
  150. FILE *fp_out = optional output file
  151.  
  152. on return from <read()>:
  153. zero nquad_skip 
  154. set quad_cycle = quad_test for use if quad_test >1
  155.  
  156. in loop calling <simpfit()> and <simpdev()>:
  157. reset maxiter to iter+prt_cycle
  158. test and reset as appropriate nquad_skip and quad_test
  159. as appropriate execute options altering flow or test values
  160.  
  161. on return from <simpfit()>:
  162. the structure <p[nvert]> has the current simplex and <pcent> the centroid,
  163. the variables mean_func, rms_func, test, and rms_data contain the mnemonically 
  164. indicated information.
  165.  
  166.  
  167. on return from <simpdev()>:
  168. the array <qmat[nfree][nfree]> contains the variance-covariance matrix, 
  169. the elements of structures <q> and <pmin> and variables yzero, ymin, ypmin,
  170. and mse contain the mnemonically indicated information.
  171. */
  172.  
  173.  
  174.  
  175. main(argc, argv)
  176. int argc ;
  177. char **argv ;
  178. {
  179.     int j, c, itemp ;
  180.     int nquad_skip, quad_cycle;
  181.     char str_buf[20] ;    
  182.  
  183.     FIL┼ *fp_ou⌠ ;
  184.  
  185.     int read_data() ;    
  186.     int simpfit() ;
  187.     int simpdev() ;
  188.     void ffitprint(), fdatprint(), fquadprint() ;
  189.  
  190.     void exit(), printf() ;
  191.     int fclose(), getchar() ;
  192.     int strlen() ;
  193.     double atof() ;
  194.     char *gets() ;
  195.  
  196.     extern char    title[] ;
  197.     extern int     nvert, nfree, nparm, iter, maxiter, maxquad_skip ;
  198.     extern int    prt_cycle, ndatval ;
  199.     extern double     quad_test, test, ypmin, yzero ;
  200.     extern struct pstruct    p[], pcent ;
  201.  
  202.     static char *paus_mess = 
  203. "\ntype ^C to stop, ^X to re-enter simpfit, any other key " ;
  204.  
  205.                 /* read data from file
  206.                 set up io */
  207.     read_data(argc, argv, stdout, &fp_out) ;
  208.  
  209.     quad_cycle = quad_test ;
  210.     nquad_skip = 0 ; 
  211.  
  212.  
  213.  
  214.                 /* BEGIN LOOP */
  215.                 /* simplex minimization alternates with 
  216.                 quadratic fit, until select exit */
  217.     while (TRUE) {
  218.         maxiter = iter + prt_cycle ;
  219.  
  220.  
  221.  
  222.                 /* carry out simplex minimization */
  223.         simpfit(stdout) ;
  224.  
  225.  
  226.                 /* print summary of simplex fitting 
  227.                 after pagσ ejec⌠ */
  228.         ffitprint(stdout) ;
  229.         if (fp_out != NULL)
  230.             ffitprint(fp_out) ;
  231.  
  232.         
  233.                 /* if test vs quad_test false, 
  234.                 loop back to simpfit ;
  235.                 the iter vs maxiter test is to allow bypass 
  236.                 on KBHIT (operator) exit from simpfit */
  237.         if (iter == maxiter) {
  238.             if (quad_test >= 1) {
  239.                 if (iter < quad_test)
  240.                     continue ;
  241.             } else if (test > quad_test)
  242.                     continue ;
  243.         } 
  244.  
  245.  
  246.                 /* clear keyboard */
  247.         while (KBHIT) {
  248.             getchar() ;
  249.         }
  250.  
  251.  
  252.                 /* if operator exit from simpfit()
  253.                 pause before continue */
  254.                 /* option to fix parameter, */
  255.                 /* or to exit, or to return to fitting */
  256. loop1:
  257.         if (iter != maxiter) {
  258.             printf("%sto display data\n", paus_mess) ;
  259.             printf("     also ^F to fix a parameter\n") ;
  260.             if ((c = getchar()) == CTRLC)
  261.                 break ;
  262.             else if (c == CTRLX)
  263.                 continue ;
  264.             else if (c == CTRLF) {
  265.                 printf("\nenter parameter (0 to %d) to fix: ",
  266.                      (nparm - 1)) ;
  267.                 if ((c = (getchar() - 0x30)) < 0 || c >= nparm)
  268.                     goto loop1 ;
  269.                 itemp = 0 ;
  270.                 for (j = 0; j < nvert; j++)
  271.                         iµ (ABS(p[j].parm[c▌ - pcent.parm[c]⌐ 
  272.                                 < 1.e-16)
  273.                         itemp++ ;
  274.                 if (itemp != nvert) {
  275.                     nfree = nfree - 1 ;
  276.                     nvert = nvert - 1 ;
  277.                 }
  278.                 else printf(
  279. "\nparameter(%d) already fixed at %17.11e", c, pcent.parm[c]) ;
  280.                 printf(
  281. "\nfix at entered value or <CR> => %17.11e ", pcent.parm[c]) ;
  282.                 if (strlen(gets(str_buf)) > 0)
  283.                     pcent.parm[c] = atof(str_buf) ;
  284.                 for (j = 0; j < nvert; j++)
  285.                     p[j].parm[c] = pcent.parm[c] ;
  286.                 goto loop1 ;
  287.             }
  288.         }
  289.     
  290.  
  291.                 /* print data array 
  292.                 after page eject */
  293.  
  294.         fdatprint(stdout) ;
  295.         if (fp_out != NULL)
  296.             fdatprint(fp_out) ;
  297.  
  298.     
  299.                 /* if operator exit from simpfit()
  300.                 pause before continue */
  301.         if (iter != maxiter) {
  302.             printf("%sfor quadratic fit\n", paus_mess) ;
  303.             if ((c = getchar()) == CTRLC)
  304.                 break ;
  305.             else if (c == CTRLX)
  306.                 continue ;
  307.         }
  308.     
  309.  
  310.                 /* carry out quadratic fit */
  311.         simpdev(stdout) ;
  312.  
  313.  
  314.                 /* print summary of results of quad fit 
  315.                 after page eject */
  316.         fquadprint(stdout) ;
  317.         if (fp_out != NULL)
  318.             fquadprint(fp_out) ;
  319.  
  320.  
  321.                 /* increment nquad_skip if ypmin > yzero */
  322.         if (ypmin < yzero) 
  323.             nquad_skip = 0 ;
  324.         else if (nquad_skip < maxquad_skip)
  325.             nquad_skip++ ;
  326.  
  327.  
  328.                 /* alter quad_test for next set of 
  329.                 fitting cycles according to nquad_skip = the
  330.                 number of quadratic fit failures and quad_test
  331.                 greater or less than unity*/
  332.         if (quad_test >= 1) {
  333.             if (iter >= quad_test)
  334.                 quad_test = iter + 
  335.                     (nquad_skip + 1) * quad_cycle ;
  336.         } else if (test <= quad_test)
  337.                 quad_test = quad_test / 10 ;
  338.  
  339.  
  340.                 /* if operator exit from simpfit()
  341.                 pause before continue */
  342.         if (iter != maxiter) {
  343.             printf("%salso re-enters simpfit\n", paus_mess) ;
  344.             if ((c = getchar()) == CTRLC)
  345.                 break ;
  346.         }
  347.     }            /* END OF LOOP */
  348.  
  349.  
  350. if (fp_out != NULL)
  351.         fclose(fp_out) ;
  352.     exit (OK) ;
  353. }                /* END OF MAIN                */
  354.