src/gmxlib/xvgr.c

   1 /*
   2  * $Id$
   3  *
   4  *       This source code is part of
   5  *
   6  *        G   R   O   M   A   C   S
   7  *
   8  * GROningen MAchine for Chemical Simulations
   9  *
  10  *               VERSION 2.0
  11  *
  12  * Copyright (c) 1991-1999
  13  * BIOSON Research Institute, Dept. of Biophysical Chemistry
  14  * University of Groningen, The Netherlands
  15  *
  16  * Please refer to:
  17  * GROMACS: A message-passing parallel molecular dynamics implementation
  18  * H.J.C. Berendsen, D. van der Spoel and R. van Drunen
  19  * Comp. Phys. Comm. 91, 43-56 (1995)
  20  *
  21  * Also check out our WWW page:
  22  * http://md.chem.rug.nl/~gmx
  23  * or e-mail to:
  24  * gromacs@chem.rug.nl
  25  *
  26  * And Hey:
  27  * Green Red Orange Magenta Azure Cyan Skyblue
  28  */
  29 static char *SRCID_xvgr_c = "$Id$";
  30
  31 #include <string.h>
  32 #include <ctype.h>
  33 #include "sysstuff.h"
  34 #include "string2.h"
  35 #include "futil.h"
  36 #include "statutil.h"
  37 #include "copyrite.h"
  38 #include "smalloc.h"
  39 #include "xvgr.h"
  40
  41 bool bXmGrace(void)
  42 {
  43   return (getenv("XMGRACE") != NULL);
  44 }
  45
  46 FILE *xgopen(char *fn,char *header,char *title,char *xaxis,char *yaxis)
  47 {
  48   FILE *xg;
  49
  50   xg=ffopen(fn,"w");
  51   fprintf(xg,"TitleText: %s: %s\n",header,title);
  52   fprintf(xg,"XUnitText: %s\n",xaxis);
  53   fprintf(xg,"YUnitText: %s\n",yaxis);
  54
  55   return xg;
  56 }
  57
  58 FILE *xvgropen(char *fn,char *title,char *xaxis,char *yaxis)
  59 {
  60   FILE *xvgr;
  61   time_t t;
  62
  63   xvgr=(FILE *)ffopen(fn,"w");
  64   fprintf(xvgr,"# This file was created by %s\n",Program());
  65   fprintf(xvgr,"# which is part of G R O M A C S:\n");
  66   fprintf(xvgr,"# %s\n",bromacs());
  67   time(&t);
  68   fprintf(xvgr,"# All this happened at: %s",ctime(&t));
  69   fprintf(xvgr,"#\n");
  70   fprintf(xvgr,"@    title \"%s\"\n",title);
  71   fprintf(xvgr,"@    xaxis  label \"%s\"\n",xaxis);
  72   fprintf(xvgr,"@    yaxis  label \"%s\"\n",yaxis);
  73   fprintf(xvgr,"@TYPE nxy\n");
  74
  75   return xvgr;
  76 }
  77
  78 void xgraph_file(char *fn,char *opts)
  79 {
  80   char doit[STRLEN];
  81
  82   if (bDoView()) {
  83     if (!opts)
  84       opts="";
  85     sprintf(doit,"xgraph =100+100 %s %s &",opts,fn);
  86     system(doit);
  87   }
  88 }
  89
  90 void xvgr_file(char *fn,char *opts)
  91 {
  92   char   doit[STRLEN];
  93
  94   if (bDoView()) {
  95     sprintf(doit,"%s -geometry +100+100 %s %s &",
  96             bXmGrace() ? "xmgrace" : "xmgr",opts ? opts : "",fn);
  97     system(doit);
  98   }
  99 }
 100
 101 void xv_file(char *fn,char *opts)
 102 {
 103   char   doit[STRLEN];
 104
 105   if (bDoView()) {
 106     sprintf(doit,"xv %s %s &",opts ? opts : "",fn);
 107     system(doit);
 108   }
 109 }
 110
 111 void xvgr_view(FILE *out,real xmin,real ymin,real xmax,real ymax)
 112 {
 113   fprintf(out,"@ view %g, %g, %g, %g\n",xmin,ymin,xmax,ymax);
 114 }
 115
 116 void xvgr_world(FILE *out,real xmin,real ymin,real xmax,real ymax)
 117 {
 118   fprintf(out,"@ world xmin %g\n"
 119           "@ world ymin %g\n"
 120           "@ world xmax %g\n"
 121           "@ world ymax %g\n",xmin,ymin,xmax,ymax);
 122 }
 123
 124 void xvgr_legend(FILE *out,int nsets,char *setname[])
 125 {
 126   int i;
 127
 128   xvgr_view(out,0.15,0.15,0.75,0.85);
 129   fprintf(out,"@ legend on\n");
 130   fprintf(out,"@ legend box on\n");
 131   fprintf(out,"@ legend loctype view\n");
 132   fprintf(out,"@ legend %g, %g\n",0.78,0.8);
 133   fprintf(out,"@ legend length %d\n",2);
 134   for(i=0; (i<nsets); i++)
 135     if (bXmGrace())
 136       fprintf(out,"@ s%d legend \"%s\"\n",i,setname[i]);
 137     else
 138       fprintf(out,"@ legend string %d \"%s\"\n",i,setname[i]);
 139 }
 140
 141 void xvgr_line_props(FILE *out, int NrSet, int LineStyle, int LineColor)
 142 {
 143   fprintf(out, "@    with g0\n");
 144   fprintf(out, "@    s%d linestyle %d\n", NrSet, LineStyle);
 145   fprintf(out, "@    s%d color %d\n", NrSet, LineColor);
 146 }
 147
 148 static char *LocTypeStr[] = { "view", "world" };
 149 static char *BoxFillStr[] = { "none", "color", "pattern" };
 150
 151 void xvgr_box(FILE *out,
 152               int LocType,
 153               real xmin,real ymin,real xmax,real ymax,
 154               int LineStyle,int LineWidth,int LineColor,
 155               int BoxFill,int BoxColor,int BoxPattern)
 156 {
 157   fprintf(out,"@with box\n");
 158   fprintf(out,"@    box on\n");
 159   fprintf(out,"@    box loctype %s\n",LocTypeStr[LocType]);
 160   fprintf(out,"@    box %g, %g, %g, %g\n",xmin,ymin,xmax,ymax);
 161   fprintf(out,"@    box linestyle %d\n",LineStyle);
 162   fprintf(out,"@    box linewidth %d\n",LineWidth);
 163   fprintf(out,"@    box color %d\n",LineColor);
 164   fprintf(out,"@    box fill %s\n",BoxFillStr[BoxFill]);
 165   fprintf(out,"@    box fill color %d\n",BoxColor);
 166   fprintf(out,"@    box fill pattern %d\n",BoxPattern);
 167   fprintf(out,"@box def\n");
 168 }
 169
 170 void lsq_y_ax(int n, real x[], real y[], real *a)
 171 {
 172   int    i;
 173   double xx,yx;
 174
 175   yx=xx=0.0;
 176   for (i=0; (i < n); i++) {
 177     yx+=y[i]*x[i];
 178     xx+=x[i]*x[i];
 179   }
 180   *a=yx/xx;
 181 }
 182
 183 void lsq_y_ax_b(int n, real x[], real y[], real *a, real *b)
 184 {
 185   int    i;
 186   double yx,xx,sx,sy;
 187
 188   yx=xx=sx=sy=0.0;
 189   for (i=0; (i < n); i++) {
 190     yx+=y[i]*x[i];
 191     xx+=x[i]*x[i];
 192     sx+=x[i];
 193     sy+=y[i];
 194   }
 195   *a=(n*yx-sy*sx)/(n*xx-sx*sx);
 196   *b=(sy-(*a)*sx)/n;
 197 }
 198
 199 static char *fgets3(FILE *fp)
 200 {
 201   static char *ptr = NULL;
 202   static int  len  = STRLEN;
 203   char *p;
 204   int  slen;
 205
 206   if (ptr==NULL)
 207     snew(ptr,len);
 208
 209   if (fgets(ptr,len-1,fp) == NULL)
 210     return NULL;
 211   p = ptr;
 212   while ((strchr(ptr,'\n') == NULL) && (!feof(fp))) {
 213     /* This line is longer than len characters, let's increase len! */
 214     len += STRLEN;
 215     p   += STRLEN;
 216     srenew(ptr,len);
 217     if (fgets(p-1,STRLEN,fp) == NULL)
 218       break;
 219   }
 220   slen = strlen(ptr);
 221   if (ptr[slen-1] == '\n')
 222     ptr[slen-1] = '\0';
 223   return ptr;
 224 }
 225
 226 static int wordcount(char *ptr)
 227 {
 228   int i,n,is[2];
 229   int cur=0;
 230 #define prev (1-cur)
 231
 232   if (strlen(ptr) == 0)
 233     return 0;
 234   /* fprintf(stderr,"ptr='%s'\n",ptr); */
 235   n=1;
 236   for(i=0; (ptr[i] != '\0'); i++) {
 237     is[cur] = isspace(ptr[i]);
 238     if ((i > 0)  && (is[cur] && !is[prev]))
 239       n++;
 240     cur=prev;
 241   }
 242   return n;
 243 }
 244
 245 int read_xvg(char *fn,real ***y,int *ny)
 246 {
 247   FILE   *fp;
 248   char   *ptr;
 249   char   *base=NULL;
 250   char   *fmt=NULL;
 251   int    k,line=0,nny,nx,maxx,rval;
 252   double lf;
 253   real   **yy=NULL;
 254
 255   *ny  = 0;
 256   nny  = 0;
 257   nx   = 0;
 258   maxx = 0;
 259   fp   = ffopen(fn,"r");
 260   while ((ptr = fgets3(fp)) != NULL) {
 261     line++;
 262     trim(ptr);
 263     if ((ptr[0] != '@') && (ptr[0] != '#')) {
 264       if (nny == 0) {
 265         (*ny) = nny = wordcount(ptr);
 266         /* fprintf(stderr,"There are %d columns in your file\n",nny);*/
 267         if (nny == 0)
 268           return 0;
 269         snew(yy,nny);
 270         snew(fmt,3*nny+1);
 271         snew(base,3*nny+1);
 272       }
 273       /* Allocate column space */
 274       if (nx >= maxx) {
 275         maxx+=1024;
 276         for(k=0; (k<nny); k++)
 277           srenew(yy[k],maxx);
 278       }
 279       /* Initiate format string */
 280       fmt[0]  = '\0';
 281       base[0] = '\0';
 282
 283       /* fprintf(stderr,"ptr='%s'\n",ptr);*/
 284       for(k=0; (k<nny); k++) {
 285         strcpy(fmt,base);
 286         strcat(fmt,"%lf");
 287         rval = sscanf(ptr,fmt,&lf);
 288         /* fprintf(stderr,"rval = %d\n",rval);*/
 289         if ((rval == EOF) || (rval == 0))
 290           break;
 291         yy[k][nx] = lf;
 292         srenew(fmt,3*(nny+1)+1);
 293         srenew(base,3*nny+1);
 294         strcat(base,"%*s");
 295       }
 296       if (k != nny) {
 297         fprintf(stderr,"Only %d columns on line %d in file %s\n",
 298                 k,line,fn);
 299         for( ; (k<nny); k++)
 300           yy[k][nx] = 0.0;
 301       }
 302       nx++;
 303     }
 304   }
 305   ffclose(fp);
 306
 307   *y = yy;
 308
 309   return nx;
 310 }
 311
 312 void dump_xvg(char *fn,char *title,int nx,int ny,real **y)
 313 {
 314   FILE *fp;
 315   int  i,j;
 316
 317   fp=xvgropen(fn,title,"X","Y");
 318   for(i=0; (i<nx); i++) {
 319     for(j=0; (j<ny); j++) {
 320       fprintf(fp,"  %12.5e",y[j][i]);
 321     }
 322     fprintf(fp,"\n");
 323   }
 324   fclose(fp);
 325 }
 326