src/tools/sas2mat.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  * Great Red Oystrich Makes All Chemists Sane
  28  */
  29 static char *SRCID_sas2mat_c = "$Id$";
  30
  31 #include <string.h>
  32 #include "sysstuff.h"
  33 #include "matio.h"
  34 #include "copyrite.h"
  35 #include "macros.h"
  36 #include "statutil.h"
  37 #include "smalloc.h"
  38
  39 int main(int argc,char *argv[])
  40 {
  41   static char *desc[] = {
  42     "sas2mat converts matrix data in [IT]raw[it] format to X PixMap format,",
  43     "which can be digested by xpm2ps to make nice plots.",
  44     "These [IT]raw[it] data may be generated by g_rms, do_dssp or your",
  45     "own program.[PAR]",
  46     "The program prompts the user for some parameters:[PAR]",
  47     "[TT]Enter nres, res0, nframes, dt, t0, nlevels:[tt][PAR]",
  48     "In this context nres is the number of residues, res0 the starting residue",
  49     "dt is the time step, t0 is the starting time, nlevels is the number",
  50     "of levels for coloring. By default a greyscale colormap is generated."
  51   };
  52   static bool   bCol=FALSE;
  53   static char   *title="Area (nm^2)";
  54   static real   ssmin=-1,ssmax=-1,t0=0,dt=1;
  55   static int    nres=1,nframes=1,r0=0,nlevels=20,nskip=0;
  56   t_pargs pa[] = {
  57     { "-col",     FALSE,  etBOOL, &bCol,
  58       "The user is prompted for rgb lower and upper values" },
  59     { "-min",     FALSE,  etREAL, &ssmin,
  60       "Lower values for the data, calculated from the data by default" },
  61     { "-max",     FALSE,  etREAL, &ssmax,
  62       "Upper values for the data, see above" },
  63     { "-title",   FALSE,  etSTR,  &title,
  64       "Title for the graph" },
  65     { "-nlevel",  FALSE,  etINT,  &nlevels,
  66       "Number of levels in graph" },
  67     { "-nres",    FALSE,  etINT,  &nres,
  68       "Number of residues (Y-axis)" },
  69     { "-nframes", FALSE,  etINT,  &nframes,
  70       "Number of frames (Y-axis)" },
  71     { "-res0",    FALSE,  etINT,  &r0,
  72       "Number of first residue" },
  73     { "-nskip",   FALSE,  etINT,  &nskip,
  74       "Number of frames to skip after every frame" },
  75     { "-dt",      FALSE,  etREAL, &dt,
  76       "Time between time frames" },
  77     { "-t0",      FALSE,  etREAL, &t0,
  78       "Time of first time frame" }
  79   };
  80
  81   FILE   *in,*out;
  82   int    i,j,k,ihi;
  83   double s;
  84   real   **ss,lo,hi,s1min,s1max;
  85   real   *resnr,*t;
  86   bool   bCheck=TRUE;
  87   t_rgb  rlo,rhi;
  88   t_filenm fnm[] = {
  89     { efOUT, "-f", "area", ffREAD },
  90     { efXPM, "-o", "sas",  ffWRITE }
  91   };
  92 #define NFILE asize(fnm)
  93
  94   /* If we want to read all frames nskip must be greater than zero */
  95   nskip += 1;
  96
  97   CopyRight(stdout,argv[0]);
  98
  99   parse_common_args(&argc,argv,0,TRUE,NFILE,fnm,asize(pa),pa,asize(desc),desc,
 100                     0,NULL);
 101
 102   snew(ss,nres);
 103   snew(resnr,nres);
 104   snew(t,nframes);
 105   for(i=0; (i<nframes); i++)
 106     t[i]=t0+i*dt;
 107   for(i=0; (i<nres); i++) {
 108     snew(ss[i],nframes);
 109   }
 110   in=ftp2FILE(efOUT,NFILE,fnm,"r");
 111   for(i=k=0; (i<nframes); i++) {
 112     for(j=0; (j<nres); j++) {
 113       fscanf(in,"%lf",&s);
 114       ss[j][k]=s;
 115     }
 116     if (!nskip || ((i % nskip) == 0))
 117       k++;
 118   }
 119   fclose(in);
 120   nframes=k;
 121
 122   lo=10000;
 123   hi=0;
 124   for(j=0; (j<nres); j++) {
 125     /* Find lowest SAS value and subtract that from all occurrences */
 126     s1min=10000;
 127     s1max=0;
 128     for(i=0; (i<nframes); i++) {
 129       s1min=min(s1min,ss[j][i]);
 130       s1max=max(s1max,ss[j][i]);
 131     }
 132     printf("res %d: ssmin=%g, ssmax=%g, diff=%g\n",j,s1min,s1max,s1max-s1min);
 133     hi=max(hi,s1max);
 134     lo=min(lo,s1min);
 135   }
 136   printf("Lowest and Highest SAS value: %g %g\n",lo,hi);
 137
 138   if (ssmin == -1)
 139     ssmin=lo;
 140   if (ssmax == -1)
 141     ssmax=hi;
 142
 143   /*
 144     hi=ssmax-ssmin;
 145     for(j=0; (j<nres); j++) {
 146     for(i=0; (i<nframes); i++)
 147     ss[j][i]-=ssmin;
 148     }
 149     */
 150
 151   /* ihi=hi; */
 152   rhi.r=0,rhi.g=0,rhi.b=0;
 153   rlo.r=1,rlo.g=1,rlo.b=1;
 154   if (bCol) {
 155     printf("Color entries:\n""drlo glo blo rhi ghi bhi\n");
 156     scanf("%f%f%f%f%f%f",&rlo.r,&rlo.g,&rlo.b,&rhi.r,&rhi.g,&rhi.b);
 157   }
 158   /*
 159   write_mapfile(ftp2fn(efMAP,NFILE,fnm),&nlevels,rlo,rhi,ssmin,ssmax);
 160   */
 161
 162   for(i=0;i<nres;i++)
 163     resnr[i]=i+1;
 164   out=ftp2FILE(efXPM,NFILE,fnm,"w");
 165   /*
 166   write_matrix(out,nres,nframes,resnr,t,ss,NULL,title,0,hi,nlevels);
 167   */
 168   write_xpm(out,"????","????","Time (ps)","Residue",
 169             nres,nframes,resnr,t,ss,ssmin,ssmax,rlo,rhi,&nlevels);
 170
 171   thanx(stdout);
 172
 173   return 0;
 174 }