home *** CD-ROM | disk | FTP | other *** search
/ OS/2 Shareware BBS: Science / Science.zip / lsqrft15.zip / plot.c < prev    next >
Text File  |  1994-01-01  |  10KB  |  389 lines
  1. #include <stdio.h>
  2. #include <stdlib.h>
  3. #include "fit.h"
  4.  
  5. extern FILE *pstream;
  6. extern double *xmin, *xmax;
  7. extern int wiflag;
  8. extern int debug;
  9. extern int noset;
  10.  
  11. int plot(func,data, order, num_indep, ndata, filename, comment, a, ma)
  12. int (*func)();
  13. double **data;
  14. struct data_order order;
  15. int num_indep;
  16. int ndata;
  17. char *filename;
  18. char *comment;
  19. double *a;
  20. int ma;
  21. {
  22. int i, failed;
  23. FILE *stream;
  24. char gnubuf[100];
  25. double umin = 1e200, umax= -1e200, vmin = 1e200, vmax=-1e200;
  26. double chisqr;
  27. int num_samples;
  28. char fitfile[40];
  29.  
  30. #ifdef DOS
  31. gnucmd(gnubuf);
  32. return 1;
  33. #endif
  34. #ifdef UNIX
  35. strcpy(fitfile,"/tmp/fit.tmp");
  36. #endif
  37. #ifndef UNIX
  38. strcpy(fitfile,"fit.tmp");
  39. #endif
  40.  
  41. if(ndata < 1 ){
  42.     printf("pl failed, no data\n");
  43.     return 0;
  44. }
  45. if( num_indep > 2){
  46.     printf("pl failed, can only plot 1 or 2 independent variables\n");
  47.     return 0;
  48. }
  49.  
  50. sprintf(gnubuf,"clear\n");
  51. gnucmd( gnubuf);
  52.  
  53. if(num_indep == 1 && !noset){
  54.     sprintf(gnubuf,"set aut\n");
  55.     gnucmd( gnubuf);
  56.     sprintf(gnubuf,"set xlabel 'x'\n");
  57.     gnucmd( gnubuf);
  58.     if(wiflag){
  59.         sprintf(gnubuf,"set xrange [%g:%g]\n", xmin[0], xmax[0]);
  60.         gnucmd( gnubuf);
  61.     }
  62.     sprintf(gnubuf,"set ylabel 'f(x)'\n");
  63.     gnucmd( gnubuf);
  64.     sprintf(gnubuf,"set title 'data and fit'\n");
  65.     gnucmd( gnubuf);
  66.     }
  67.  
  68.     if(num_indep == 2 && !noset){
  69.     sprintf(gnubuf,"set aut\n");
  70.     gnucmd( gnubuf);
  71.     sprintf(gnubuf,"set xlabel 'x'\n");
  72.     gnucmd( gnubuf);
  73.     if(wiflag){
  74.         sprintf(gnubuf,"set xrange [%g:%g]\n", xmin[0], xmax[0]);
  75.         gnucmd( gnubuf);
  76.         sprintf(gnubuf,"set yrange [%g:%g]\n", xmin[1], xmax[1]);
  77.         gnucmd( gnubuf);
  78.     }
  79.     sprintf(gnubuf,"set ylabel 'y'\n");
  80.     gnucmd( gnubuf);
  81.     sprintf(gnubuf,"set zlabel 'f(x,y)'\n");
  82.     gnucmd( gnubuf);
  83.     sprintf(gnubuf,"set title 'data and fit'\n");
  84.     gnucmd( gnubuf);
  85. }
  86.  
  87.  
  88. /* function is plotted through a defined gnuplot function */
  89. if(strncmp(comment,"FILE",4) != 0){
  90.  
  91.         /* define the parameters in gnuplot */
  92.     for(i = 0; i < ma; i++){
  93.         sprintf(gnubuf,"a%d = %g\n",i,a[i]);
  94.         gnucmd(gnubuf);
  95.     }
  96.  
  97.         /* define the function in gnuplot */
  98.     sprintf(gnubuf,"%s\n",comment);
  99.     gnucmd( gnubuf);
  100.  
  101.         /* tell gnuplot to do the plot */
  102.     if(num_indep == 1){
  103.         sprintf(gnubuf,"plot '%s' using %d:%d, f(x) \n",
  104.             filename,order.x[0] + 1, order.y + 1);
  105.         gnucmd( gnubuf);
  106.         if(!noset){
  107.             if(ndata > 40) num_samples = ndata/2;
  108.             else if(ndata > 20) num_samples = ndata;
  109.             else num_samples = 20;
  110.             if(num_samples > 500) num_samples = 500;
  111.             sprintf(gnubuf,"set samples %d\n",num_samples);
  112.             gnucmd(gnubuf);
  113.         }
  114.     }
  115.  
  116. /* for 2 independent variables, we need to use the gnuplot */
  117. /* command splot and plot parametrically.  See gnuplot documentation */
  118. /* for the reasons why.  We also need to know tell gnuplot */
  119. /* the ranges of the independent variables, so we calculate these also */
  120.  
  121.     if(num_indep == 2){
  122.         if(wiflag == 0)
  123.             for(i = 0; i < ndata; i++){
  124.                 if(data[order.x[0]][i] < umin) umin = data[order.x[0]][i];
  125.                 if(data[order.x[0]][i] > umax) umax = data[order.x[0]][i];
  126.                 if(data[order.x[1]][i] < vmin) vmin = data[order.x[1]][i];
  127.                 if(data[order.x[1]][i] > vmax) vmax = data[order.x[1]][i];
  128.             }
  129.         else{
  130.             umin = xmin[0];
  131.             umax = xmax[0];
  132.             vmin = xmin[1];
  133.             vmax = xmax[1];
  134.         }
  135.         sprintf(gnubuf,"set parametric\n");
  136.         gnucmd( gnubuf);
  137.         if(1){
  138.             sprintf(gnubuf,"set urange [%g:%g]\n", umin, umax);
  139.             gnucmd( gnubuf);
  140.             sprintf(gnubuf,"set vrange [%g:%g]\n", vmin, vmax);
  141.             gnucmd(gnubuf);
  142.             sprintf(gnubuf,"set samples 40\n");
  143.             gnucmd(gnubuf);
  144.         }
  145.         sprintf(gnubuf,"splot '%s' using %d:%d:%d, u,v,f(u,v) \n",
  146.             filename,order.x[0] + 1, order.x[1] + 1, order.y + 1);
  147.         gnucmd(gnubuf);
  148.         sprintf(gnubuf,"set noparametric\n");
  149.         gnucmd(gnubuf);
  150.     }
  151.  
  152. }
  153.  
  154. /* function is plotted through a temporary file */
  155. else if(strncmp(comment,"FILE",4) == 0){
  156.     if(debug) printf("in plot(), about to call calc_yfit()\n");
  157.  
  158. /* calculate the values of the fitting function for current parameters */
  159.     failed = calc_yfit(func, data, order, num_indep, a, ndata, ma, &chisqr);
  160.     if(failed) return 1;
  161.     if(debug) printf("in plot(), returned from calc_yfit()\n");
  162.     stream = fopen(fitfile,"w");
  163.     if(num_indep == 1)
  164.         for(i = 0; i < ndata; i++) 
  165.             if(wiflag == 0 || (data[order.x[0]][i] >= xmin[0] && data[order.x[0]][i] <= xmax[0]))
  166.                 fprintf(stream,"%g %g\n", data[order.x[0]][i], data[order.yfit][i]);
  167.     if(num_indep == 2)
  168.         for(i = 0; i < ndata; i++)
  169.             if(wiflag == 0 || (data[order.x[0]][i] >= xmin[0] && data[order.x[0]][i] <= xmax[0]) && (data[order.x[1]][i] >= xmin[1] && data[order.x[1]][i] <= xmax[1]))
  170.                 fprintf(stream,"%g %g %g\n", data[order.x[0]][i], data[order.x[1]][i], data[order.yfit][i]);
  171.     fclose(stream);
  172.     if(num_indep == 1){
  173.         sprintf(gnubuf,"plot '%s' using %d:%d,'%s'\n",
  174.                                 filename,order.x[0] + 1, order.y + 1, fitfile);
  175.         gnucmd(gnubuf);
  176.     }
  177.     if(num_indep == 2){
  178.         sprintf(gnubuf,"set parametric\n");
  179.         gnucmd(gnubuf);
  180.         sprintf(gnubuf,"splot '%s' using %d:%d:%d,'%s' \n",
  181.                                 filename,order.x[0] + 1, order.x[1] + 1,order.y + 1, fitfile);
  182.         gnucmd(gnubuf);
  183.         sprintf(gnubuf,"set noparametric\n");
  184.         gnucmd(gnubuf);
  185.     }
  186.  
  187. }
  188. return 0;
  189. }
  190.  
  191. /* the pr function plot residual errors */
  192. int pr(command,func,data, order, num_indep, ndata, filename, comment, a, ma)
  193. char *command;
  194. int (*func)();
  195. double **data;
  196. struct data_order order;
  197. int num_indep;
  198. int ndata;
  199. char *filename;
  200. char *comment;
  201. double *a;
  202. int ma;
  203. {
  204. int i, flag, failed;
  205. FILE *stream;
  206. char gnubuf[100];
  207. double chisqr;
  208. char fitfile[40];
  209.  
  210. #ifdef DOS
  211. gnucmd(gnubuf);
  212. return 1;
  213. #endif
  214. #ifdef UNIX
  215. strcpy(fitfile,"/tmp/fit.tmp");
  216. #endif
  217. #ifndef UNIX
  218. strcpy(fitfile,"fit.tmp");
  219. #endif
  220.  
  221. if(debug)printf("in pr() command: %s\n", command);
  222. if(ndata < 1 ){
  223.     printf("pr failed, no data\n");
  224.     return 0;
  225. }
  226. if( num_indep > 2){
  227.     printf("pr failed, can only plot 1 or 2 independent variables\n");
  228.     return 0;
  229. }
  230.  
  231. sprintf(gnubuf,"clear\n");
  232. if( sscanf(command,"%*s %d", &flag ) == 0) flag = 1;
  233. if(debug)printf("flag: %d\n", flag);
  234.  
  235. failed = calc_yfit(func, data, order, num_indep, a, ndata, ma, &chisqr);
  236. if(failed) return 1;
  237.  
  238. /* option 1, plot yfit vs. y */
  239. if(flag ==1){
  240.     stream = fopen(fitfile,"w");
  241.     for(i = 0; i < ndata; i++) 
  242.         if(!wiflag || (data[order.x[0]][i] >= xmin[0] && data[order.x[0]][i] <= xmax[0])) fprintf(stream,"%g %g\n", data[order.yfit][i], data[order.y][i]);
  243.     fclose(stream);
  244.     if(!noset){
  245.         sprintf(gnubuf,"set aut\n");
  246.         gnucmd( gnubuf);
  247.         sprintf(gnubuf,"set title 'data vs. fit'\n");
  248.         gnucmd( gnubuf);
  249.         sprintf(gnubuf,"set xlabel 'f(x)'\n");
  250.         gnucmd( gnubuf);
  251.         if(wiflag){
  252.             sprintf(gnubuf,"set xrange [%g:%g]\n", xmin[0], xmax[0]);
  253.             gnucmd( gnubuf);
  254.         }
  255.         sprintf(gnubuf,"set ylabel 'y'\n");
  256.         gnucmd( gnubuf);
  257.     }
  258.     sprintf(gnubuf,"plot '%s' with points \n", fitfile);
  259.     gnucmd(gnubuf);
  260. }
  261.  
  262. /* option 2, plot y-yfit as a function of independent variables */
  263. else{
  264.     if(debug)printf("stream = fopen();\n");
  265.     stream = fopen(fitfile,"w");
  266.     if(debug)printf("if(num_indep ==1);\n");
  267.     if(!noset){
  268.         sprintf(gnubuf,"set aut\n");
  269.         gnucmd( gnubuf);
  270.         sprintf(gnubuf,"set title 'data-fit'\n");
  271.         gnucmd( gnubuf);
  272.     }
  273.     if(num_indep ==1){
  274.         for(i = 0; i < ndata; i++)
  275.             if(!wiflag || (data[order.x[0]][i] >= xmin[0] && data[order.x[0]][i] <= xmax[0]))
  276.                 fprintf(stream,"%g %g\n", data[order.x[0]][i],
  277.                         data[order.y][i]-data[order.yfit][i]);
  278.     }
  279.     else if(num_indep ==2){
  280.         for(i = 0; i < ndata; i++)
  281.             if(!wiflag || (data[order.x[0]][i] >= xmin[0] && data[order.x[0]][i] <= xmax[0]))
  282.                 fprintf(stream,"%g %g %g\n", data[order.x[0]][i],data[order.x[1]][i],
  283.                         data[order.y][i]-data[order.yfit][i]);
  284.     }
  285.     fclose(stream);
  286.     if(num_indep ==1){
  287.         if(!noset){
  288.             sprintf(gnubuf,"set aut\n");
  289.             gnucmd( gnubuf);
  290.             if(wiflag){
  291.                 sprintf(gnubuf,"set xrange [%g:%g]\n", xmin[0], xmax[0]);
  292.                 gnucmd( gnubuf);
  293.             }
  294.             sprintf(gnubuf,"set title 'y - f(x)'\n");
  295.             gnucmd( gnubuf);
  296.             sprintf(gnubuf,"set xlabel 'x'\n");
  297.             gnucmd( gnubuf);
  298.             sprintf(gnubuf,"set ylabel 'y - f(x)'\n");
  299.             gnucmd( gnubuf);
  300.         }
  301.         sprintf(gnubuf,"plot '%s'\n", fitfile);
  302.         gnucmd(gnubuf);
  303.     }    
  304.     if(num_indep ==2){
  305.         if(!noset){
  306.             sprintf(gnubuf,"set aut\n");
  307.             gnucmd( gnubuf);
  308.             sprintf(gnubuf,"set title 'y - f(x)'\n");
  309.             gnucmd( gnubuf);
  310.             sprintf(gnubuf,"set xlabel 'x'\n");
  311.             gnucmd( gnubuf);
  312.             sprintf(gnubuf,"set zlabel 'y - f(x)'\n");
  313.             gnucmd( gnubuf);
  314.             sprintf(gnubuf,"set ylabel 'y'\n");
  315.             gnucmd( gnubuf);
  316.         }
  317.         sprintf(gnubuf,"set parametric\n");
  318.         gnucmd(gnubuf);
  319.         sprintf(gnubuf,"splot '%s'\n", fitfile);
  320.         gnucmd(gnubuf);
  321.         sprintf(gnubuf,"set noparametric\n");
  322.         gnucmd(gnubuf);
  323.     }
  324. }
  325. return 0;
  326. }
  327.  
  328. #ifdef OS2
  329.  
  330. int myplot(func,data, order, num_indep, ndata, filename, comment, a, ma)
  331. int (*func)();
  332. double **data;
  333. struct data_order order;
  334. int num_indep;
  335. int ndata;
  336. char *filename;
  337. char *comment;
  338. double *a;
  339. int ma;
  340. {
  341. int i, failed, j;
  342. FILE *stream;
  343. char buf[100];
  344. double chisqr;
  345. char fitfile[40];
  346. double x;
  347.  
  348. strcpy(fitfile,"fit.tmp");
  349.  
  350. if(ndata < 1 ){
  351.     printf("pl failed, no data\n");
  352.     return 0;
  353. }
  354. if( num_indep > 1){
  355.     printf("pl failed, can only plot 1 indep variable\n");
  356.     return 0;
  357. }
  358.  
  359. /* function is plotted through a temporary file */
  360.  
  361. /* calculate the values of the fitting function for current parameters */
  362. /*
  363. if(debug) printf("in plot(), about to call calc_yfit()\n");
  364. failed = calc_yfit(func, data, order, num_indep, a, ndata, ma, &chisqr);
  365. if(failed) return 1;
  366. if(debug) printf("in myplot(), returned from calc_yfit()\n");
  367. */
  368.  
  369. stream = fopen(fitfile,"w");
  370. for(i = 0; i < ndata; i++){
  371.         x = data[order.x[0]][i];
  372.         if( !wiflag || (x > xmin[0])  && (x < xmax[0]))
  373.             fprintf(stream,"%g %g %g\n", 
  374.                     x, data[order.y][i], data[order.yfit][i]);
  375. }
  376. fclose(stream);
  377.  
  378. sprintf(buf,"gd fit.tmp 2\n");
  379. fitcmd(buf);
  380. sprintf(buf,"sy 0 L\n");
  381. fitcmd(buf);
  382. sprintf(buf,"pl\n");
  383. fitcmd(buf);
  384. return 0;
  385. }
  386.  
  387. #endif
  388.  
  389.