fixed g_anaeig incorrect multiple extremes, fixes #870
authorBerk Hess <hess@kth.se>
Mon, 5 Mar 2012 11:18:58 +0000 (12:18 +0100)
committerBerk Hess <hess@kth.se>
Mon, 5 Mar 2012 11:18:58 +0000 (12:18 +0100)
Change-Id: Id06f267c1da300f2a940778833ec35b85fecd3f1

src/tools/gmx_anaeig.c

index a32b5095eac24add2904a379a710e129f81be1d9..811812f17d92186942480d32f35392bd32907ef4 100644 (file)
@@ -408,7 +408,8 @@ static void project(const char *trajfile,t_topology *top,int ePBC,matrix topbox,
   int     nat,i,j,d,v,vec,nfr,nframes=0,snew_size,frame;
   t_trxstatus *out=NULL;
   t_trxstatus *status;
-  int     noutvec_extr,*imin,*imax;
+  int     noutvec_extr,imin,imax;
+  real    *pmin,*pmax;
   atom_id *all_at;
   matrix  box;
   rvec    *xread,*x;
@@ -615,29 +616,31 @@ static void project(const char *trajfile,t_topology *top,int ePBC,matrix topbox,
   }
   
   if (extremefile) {
+    snew(pmin,noutvec_extr);
+    snew(pmax,noutvec_extr);
     if (extreme==0) {
       fprintf(stderr,"%11s %17s %17s\n","eigenvector","Minimum","Maximum");
       fprintf(stderr,
-             "%11s %10s %10s %10s %10s\n","","value","time","value","time");
-      snew(imin,noutvec_extr);
-      snew(imax,noutvec_extr);
+             "%11s %10s %10s %10s %10s\n","","value","frame","value","frame");
+      imin = 0;
+      imax = 0;
       for(v=0; v<noutvec_extr; v++) {
        for(i=0; i<nframes; i++) {
-         if (inprod[v][i]<inprod[v][imin[v]])
-           imin[v]=i;
-         if (inprod[v][i]>inprod[v][imax[v]])
-           imax[v]=i;
+         if (inprod[v][i]<inprod[v][imin])
+           imin = i;
+         if (inprod[v][i]>inprod[v][imax])
+           imax = i;
        }
-       min=inprod[v][imin[v]];
-       max=inprod[v][imax[v]];
-       fprintf(stderr,"%7d     %10.6f %10.1f %10.6f %10.1f\n",
+       pmin[v] = inprod[v][imin];
+       pmax[v] = inprod[v][imax];
+       fprintf(stderr,"%7d     %10.6f %10d %10.6f %10d\n",
                eignr[outvec[v]]+1,
-               min,inprod[noutvec][imin[v]],max,inprod[noutvec][imax[v]]); 
+               pmin[v],imin,pmax[v],imax); 
       }
     }
     else {
-      min=-extreme;
-      max=+extreme;
+      pmin[0] = -extreme;
+      pmax[0] =  extreme;
     }
     /* build format string for filename: */
     strcpy(str,extremefile);/* copy filename */
@@ -661,12 +664,14 @@ static void project(const char *trajfile,t_topology *top,int ePBC,matrix topbox,
        for(i=0; i<natoms; i++)
          for(d=0; d<DIM; d++) 
            xread[index[i]][d] = 
-             (xav[i][d] + (min*(nextr-frame-1)+max*frame)/(nextr-1)
+             (xav[i][d] + (pmin[v]*(nextr-frame-1)+pmax[v]*frame)/(nextr-1)
              *eigvec[outvec[v]][i][d]/sqrtm[i]);
        write_trx(out,natoms,index,atoms,0,frame,topbox,xread,NULL,NULL);
       }
       close_trx(out);
     }
+    sfree(pmin);
+    sfree(pmax);
   }
   fprintf(stderr,"\n");
 }