1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
11 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
12 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
13 * Copyright (c) 2001-2004, The GROMACS development team,
14 * check out http://www.gromacs.org for more information.
16 * This program is free software; you can redistribute it and/or
17 * modify it under the terms of the GNU General Public License
18 * as published by the Free Software Foundation; either version 2
19 * of the License, or (at your option) any later version.
21 * If you want to redistribute modifications, please consider that
22 * scientific software is very special. Version control is crucial -
23 * bugs must be traceable. We will be happy to consider code for
24 * inclusion in the official distribution, but derived work must not
25 * be called official GROMACS. Details are found in the README & COPYING
26 * files - if they are missing, get the official version at www.gromacs.org.
28 * To help us fund GROMACS development, we humbly ask that you cite
29 * the papers on the package - you can find them in the top README file.
31 * For more info, check our website at http://www.gromacs.org
34 * Green Red Orange Magenta Azure Cyan Skyblue
63 static void low_print_data(FILE *fp,real time,rvec x[],int n,atom_id *index,
64 gmx_bool bDim[],const char *sffmt)
68 fprintf(fp," %g",time);
83 fprintf(fp,sffmt,x[ii][d]);
88 fprintf(fp,sffmt,norm(x[ii]));
94 static void average_data(rvec x[],rvec xav[],real *mass,
95 int ngrps,int isize[],atom_id **index)
102 for(g=0; g<ngrps; g++)
107 for(i=0; i<isize[g]; i++)
132 xav[g][d] = sum[d]/mtot;
137 /* mass=NULL, so these are forces: sum only (do not average) */
146 static void print_data(FILE *fp,real time,rvec x[],real *mass,gmx_bool bCom,
147 int ngrps,int isize[],atom_id **index,gmx_bool bDim[],
150 static rvec *xav=NULL;
158 average_data(x,xav,mass,ngrps,isize,index);
159 low_print_data(fp,time,xav,ngrps,NULL,bDim,sffmt);
163 low_print_data(fp,time,x,isize[0],index[0],bDim,sffmt);
167 static void write_trx_x(t_trxstatus *status,t_trxframe *fr,real *mass,gmx_bool bCom,
168 int ngrps,int isize[],atom_id **index)
170 static rvec *xav=NULL;
171 static t_atoms *atoms=NULL;
184 snew(atoms->atom,ngrps);
186 /* Note that atom and residue names will be the ones
187 * of the first atom in each group.
189 for(i=0; i<ngrps; i++)
191 atoms->atom[i] = fr->atoms->atom[index[i][0]];
192 atoms->atomname[i] = fr->atoms->atomname[index[i][0]];
195 average_data(fr->x,xav,mass,ngrps,isize,index);
197 fr_av.natoms = ngrps;
200 write_trxframe(status,&fr_av,NULL);
204 write_trxframe_indexed(status,fr,isize[0],index[0],NULL);
208 static void make_legend(FILE *fp,int ngrps,int isize,atom_id index[],
209 char **name,gmx_bool bCom,gmx_bool bMol,gmx_bool bDim[],
210 const output_env_t oenv)
213 const char *dimtxt[] = { " X", " Y", " Z", "" };
229 for(d=0; d<=DIM; d++)
236 sprintf(leg[j],"mol %d%s",index[i]+1,dimtxt[d]);
240 sprintf(leg[j],"%s%s",name[i],dimtxt[d]);
244 sprintf(leg[j],"atom %d%s",index[i]+1,dimtxt[d]);
250 xvgr_legend(fp,j,(const char**)leg,oenv);
259 static real ekrot(rvec x[],rvec v[],real mass[],int isize,atom_id index[])
261 static real **TCM=NULL,**L;
262 double tm,m0,lxx,lxy,lxz,lyy,lyz,lzz,ekrot;
286 for(i=0; i<isize; i++)
292 for(m=0; (m<DIM); m++)
294 xcm[m] += m0*x[j][m]; /* c.o.m. position */
295 vcm[m] += m0*v[j][m]; /* c.o.m. velocity */
296 acm[m] += m0*a0[m]; /* rotational velocity around c.o.m. */
300 for(m=0; (m<DIM); m++)
307 lxx=lxy=lxz=lyy=lyz=lzz=0.0;
308 for(i=0; i<isize; i++)
314 dx[m] = x[j][m] - xcm[m];
316 lxx += dx[XX]*dx[XX]*m0;
317 lxy += dx[XX]*dx[YY]*m0;
318 lxz += dx[XX]*dx[ZZ]*m0;
319 lyy += dx[YY]*dx[YY]*m0;
320 lyz += dx[YY]*dx[ZZ]*m0;
321 lzz += dx[ZZ]*dx[ZZ]*m0;
324 L[XX][XX] = lyy + lzz;
328 L[YY][YY] = lxx + lzz;
332 L[ZZ][ZZ] = lxx + lyy;
334 m_inv_gen(L,DIM,TCM);
336 /* Compute omega (hoeksnelheid) */
343 ocm[m] += TCM[m][n]*acm[n];
345 ekrot += 0.5*ocm[m]*acm[m];
351 static real ektrans(rvec v[],real mass[],int isize,atom_id index[])
359 for(i=0; i<isize; i++)
364 mvcom[d] += mass[j]*v[j][d];
369 return dnorm2(mvcom)/(mtot*2);
372 static real temp(rvec v[],real mass[],int isize,atom_id index[])
378 for(i=0; i<isize; i++)
381 ekin2 += mass[j]*norm2(v[j]);
384 return ekin2/(3*isize*BOLTZ);
387 static void remove_jump(matrix box,int natoms,rvec xp[],rvec x[])
394 hbox[d] = 0.5*box[d][d];
396 for(i=0; i<natoms; i++)
398 for(m=DIM-1; m>=0; m--)
400 while (x[i][m]-xp[i][m] <= -hbox[m])
404 x[i][d] += box[m][d];
407 while (x[i][m]-xp[i][m] > hbox[m])
411 x[i][d] -= box[m][d];
418 static void write_pdb_bfac(const char *fname,const char *xname,
419 const char *title,t_atoms *atoms,int ePBC,matrix box,
420 int isize,atom_id *index,int nfr_x,rvec *x,
422 gmx_bool bDim[],real scale_factor,
423 const output_env_t oenv)
431 if ((nfr_x == 0) || (nfr_v == 0))
433 fprintf(stderr,"No frames found for %s, will not write %s\n",
438 fprintf(stderr,"Used %d frames for %s\n",nfr_x,"coordinates");
439 fprintf(stderr,"Used %d frames for %s\n",nfr_v,title);
458 for(i=0; i<isize; i++)
460 svmul(scale,x[index[i]],x[index[i]]);
463 for(i=0; i<isize; i++)
465 svmul(scale,sum[index[i]],sum[index[i]]);
468 fp=xvgropen(xname,title,"Atom","",oenv);
469 for(i=0; i<isize; i++)
471 fprintf(fp,"%-5d %10.3f %10.3f %10.3f\n",i,
472 sum[index[i]][XX],sum[index[i]][YY],sum[index[i]][ZZ]);
477 for(i=0; i<isize; i++)
482 if (bDim[m] || bDim[DIM])
484 len2 += sqr(sum[index[i]][m]);
493 if (scale_factor != 0)
495 scale = scale_factor;
505 scale = 10.0/sqrt(max);
509 printf("Maximum %s is %g on atom %d %s, res. %s %d\n",
510 title,sqrt(max)/nfr_v,maxi+1,*(atoms->atomname[maxi]),
511 *(atoms->resinfo[atoms->atom[maxi].resind].name),
512 atoms->resinfo[atoms->atom[maxi].resind].nr);
514 if (atoms->pdbinfo == NULL)
516 snew(atoms->pdbinfo,atoms->nr);
520 for(i=0; i<isize; i++)
525 if (bDim[m] || bDim[DIM])
527 len2 += sqr(sum[index[i]][m]);
530 atoms->pdbinfo[index[i]].bfac = sqrt(len2)*scale;
535 for(i=0; i<isize; i++)
537 atoms->pdbinfo[index[i]].bfac = sum[index[i]][onedim]*scale;
540 write_sto_conf_indexed(fname,title,atoms,x,NULL,ePBC,box,isize,index);
544 static void update_histo(int gnx,atom_id index[],rvec v[],
545 int *nhisto,int **histo,real binwidth)
553 for(i=0; (i<gnx); i++)
555 vn = norm(v[index[i]]);
556 vnmax = max(vn,vnmax);
559 *nhisto = 1+(vnmax/binwidth);
560 snew(*histo,*nhisto);
562 for(i=0; (i<gnx); i++)
564 vn = norm(v[index[i]]);
569 fprintf(stderr,"Extending histogram from %d to %d\n",*nhisto,nnn);
572 for(m=*nhisto; (m<nnn); m++)
582 static void print_histo(const char *fn,int nhisto,int histo[],real binwidth,
583 const output_env_t oenv)
588 fp = xvgropen(fn,"Velocity distribution","V (nm/ps)","arbitrary units",
590 for(i=0; (i<nhisto); i++)
592 fprintf(fp,"%10.3e %10d\n",i*binwidth,histo[i]);
597 int gmx_traj(int argc,char *argv[])
599 const char *desc[] = {
600 "g_traj plots coordinates, velocities, forces and/or the box.",
601 "With [TT]-com[tt] the coordinates, velocities and forces are",
602 "calculated for the center of mass of each group.",
603 "When [TT]-mol[tt] is set, the numbers in the index file are",
604 "interpreted as molecule numbers and the same procedure as with",
605 "[TT]-com[tt] is used for each molecule.[PAR]",
606 "Option [TT]-ot[tt] plots the temperature of each group,",
607 "provided velocities are present in the trajectory file.",
608 "No corrections are made for constrained degrees of freedom!",
609 "This implies [TT]-com[tt].[PAR]",
610 "Options [TT]-ekt[tt] and [TT]-ekr[tt] plot the translational and",
611 "rotational kinetic energy of each group,",
612 "provided velocities are present in the trajectory file.",
613 "This implies [TT]-com[tt].[PAR]",
614 "Options [TT]-cv[tt] and [TT]-cf[tt] write the average velocities",
615 "and average forces as temperature factors to a pdb file with",
616 "the average coordinates. The temperature factors are scaled such",
617 "that the maximum is 10. The scaling can be changed with the option",
618 "[TT]-scale[tt]. To get the velocities or forces of one",
619 "frame set both [TT]-b[tt] and [TT]-e[tt] to the time of",
620 "desired frame. When averaging over frames you might need to use",
621 "the [TT]-nojump[tt] option to obtain the correct average coordinates.",
622 "If you select either of these option the average force and velocity",
623 "for each atom are written to an xvg file as well",
624 "(specified with [TT]-av[tt] or [TT]-af[tt]).[PAR]",
625 "Option [TT]-vd[tt] computes a velocity distribution, i.e. the",
626 "norm of the vector is plotted. In addition in the same graph",
627 "the kinetic energy distribution is given."
629 static gmx_bool bMol=FALSE,bCom=FALSE,bPBC=TRUE,bNoJump=FALSE;
630 static gmx_bool bX=TRUE,bY=TRUE,bZ=TRUE,bNorm=FALSE,bFP=FALSE;
631 static int ngroups=1;
632 static real scale=0,binwidth=1;
634 { "-com", FALSE, etBOOL, {&bCom},
635 "Plot data for the com of each group" },
636 { "-pbc", FALSE, etBOOL, {&bPBC},
637 "Make molecules whole for COM" },
638 { "-mol", FALSE, etBOOL, {&bMol},
639 "Index contains molecule numbers iso atom numbers" },
640 { "-nojump", FALSE, etBOOL, {&bNoJump},
641 "Remove jumps of atoms across the box" },
642 { "-x", FALSE, etBOOL, {&bX},
643 "Plot X-component" },
644 { "-y", FALSE, etBOOL, {&bY},
645 "Plot Y-component" },
646 { "-z", FALSE, etBOOL, {&bZ},
647 "Plot Z-component" },
648 { "-ng", FALSE, etINT, {&ngroups},
649 "Number of groups to consider" },
650 { "-len", FALSE, etBOOL, {&bNorm},
651 "Plot vector length" },
652 { "-fp", FALSE, etBOOL, {&bFP},
653 "Full precision output" },
654 { "-bin", FALSE, etREAL, {&binwidth},
655 "Binwidth for velocity histogram (nm/ps)" },
656 { "-scale", FALSE, etREAL, {&scale},
657 "Scale factor for pdb output, 0 is autoscale" }
659 FILE *outx=NULL,*outv=NULL,*outf=NULL,*outb=NULL,*outt=NULL;
660 FILE *outekt=NULL,*outekr=NULL;
667 int flags,nvhisto=0,*vhisto=NULL;
669 rvec *sumxv=NULL,*sumv=NULL,*sumxf=NULL,*sumf=NULL;
672 t_trxstatus *status_out=NULL;
673 gmx_rmpbc_t gpbc=NULL;
675 int nr_xfr,nr_vfr,nr_ffr;
678 atom_id **index0,**index;
681 gmx_bool bTop,bOX,bOXT,bOV,bOF,bOB,bOT,bEKT,bEKR,bCV,bCF;
682 gmx_bool bDim[4],bDum[4],bVD;
683 char *sffmt,sffmt6[1024];
684 const char *box_leg[6] = { "XX", "YY", "ZZ", "YX", "ZX", "ZY" };
688 { efTRX, "-f", NULL, ffREAD },
689 { efTPS, NULL, NULL, ffREAD },
690 { efNDX, NULL, NULL, ffOPTRD },
691 { efXVG, "-ox", "coord.xvg", ffOPTWR },
692 { efTRX, "-oxt","coord.xtc", ffOPTWR },
693 { efXVG, "-ov", "veloc.xvg", ffOPTWR },
694 { efXVG, "-of", "force.xvg", ffOPTWR },
695 { efXVG, "-ob", "box.xvg", ffOPTWR },
696 { efXVG, "-ot", "temp.xvg", ffOPTWR },
697 { efXVG, "-ekt","ektrans.xvg", ffOPTWR },
698 { efXVG, "-ekr","ekrot.xvg", ffOPTWR },
699 { efXVG, "-vd", "veldist.xvg", ffOPTWR },
700 { efPDB, "-cv", "veloc.pdb", ffOPTWR },
701 { efPDB, "-cf", "force.pdb", ffOPTWR },
702 { efXVG, "-av", "all_veloc.xvg", ffOPTWR },
703 { efXVG, "-af", "all_force.xvg", ffOPTWR }
705 #define NFILE asize(fnm)
707 CopyRight(stderr,argv[0]);
708 parse_common_args(&argc,argv,
709 PCA_CAN_TIME | PCA_TIME_UNIT | PCA_CAN_VIEW | PCA_BE_NICE,
710 NFILE,fnm,asize(pa),pa,asize(desc),desc,0,NULL,&oenv);
714 fprintf(stderr,"Interpreting indexfile entries as molecules.\n"
715 "Using center of mass.\n");
718 bOX = opt2bSet("-ox",NFILE,fnm);
719 bOXT = opt2bSet("-oxt",NFILE,fnm);
720 bOV = opt2bSet("-ov",NFILE,fnm);
721 bOF = opt2bSet("-of",NFILE,fnm);
722 bOB = opt2bSet("-ob",NFILE,fnm);
723 bOT = opt2bSet("-ot",NFILE,fnm);
724 bEKT = opt2bSet("-ekt",NFILE,fnm);
725 bEKR = opt2bSet("-ekr",NFILE,fnm);
726 bCV = opt2bSet("-cv",NFILE,fnm) || opt2bSet("-av",NFILE,fnm);
727 bCF = opt2bSet("-cf",NFILE,fnm) || opt2bSet("-af",NFILE,fnm);
728 bVD = opt2bSet("-vd",NFILE,fnm) || opt2parg_bSet("-bin",asize(pa),pa);
729 if (bMol || bOT || bEKT || bEKR)
741 sffmt = "\t" gmx_real_fullprecision_pfmt;
747 sprintf(sffmt6,"%s%s%s%s%s%s",sffmt,sffmt,sffmt,sffmt,sffmt,sffmt);
749 bTop = read_tps_conf(ftp2fn(efTPS,NFILE,fnm),title,&top,&ePBC,
751 bCom && (bOX || bOXT || bOV || bOT || bEKT || bEKR));
753 if ((bMol || bCV || bCF) && !bTop)
755 gmx_fatal(FARGS,"Need a run input file for option -mol, -cv or -cf");
760 indexfn = ftp2fn(efNDX,NFILE,fnm);
764 indexfn = ftp2fn_null(efNDX,NFILE,fnm);
767 if (!(bCom && !bMol))
771 snew(grpname,ngroups);
772 snew(isize0,ngroups);
773 snew(index0,ngroups);
774 get_index(&(top.atoms),indexfn,ngroups,isize0,index0,grpname);
783 for (i=0; i<ngroups; i++)
785 if (index0[0][i] < 0 || index0[0][i] >= mols->nr)
787 gmx_fatal(FARGS,"Molecule index (%d) is out of range (%d-%d)",
788 index0[0][i]+1,1,mols->nr);
790 isize[i] = atndx[index0[0][i]+1] - atndx[index0[0][i]];
791 snew(index[i],isize[i]);
792 for(j=0; j<isize[i]; j++)
794 index[i][j] = atndx[index0[0][i]] + j;
805 snew(mass,top.atoms.nr);
806 for(i=0; i<top.atoms.nr; i++)
808 mass[i] = top.atoms.atom[i].m;
819 flags = flags | TRX_READ_X;
820 outx = xvgropen(opt2fn("-ox",NFILE,fnm),
821 bCom ? "Center of mass" : "Coordinate",
822 output_env_get_xvgr_tlabel(oenv),"Coordinate (nm)",oenv);
823 make_legend(outx,ngroups,isize0[0],index0[0],grpname,bCom,bMol,bDim,oenv);
827 flags = flags | TRX_READ_X;
828 status_out = open_trx(opt2fn("-oxt",NFILE,fnm),"w");
832 flags = flags | TRX_READ_V;
833 outv = xvgropen(opt2fn("-ov",NFILE,fnm),
834 bCom ? "Center of mass velocity" : "Velocity",
835 output_env_get_xvgr_tlabel(oenv),"Velocity (nm/ps)",oenv);
836 make_legend(outv,ngroups,isize0[0],index0[0],grpname,bCom,bMol,bDim,oenv);
839 flags = flags | TRX_READ_F;
840 outf = xvgropen(opt2fn("-of",NFILE,fnm),"Force",
841 output_env_get_xvgr_tlabel(oenv),"Force (kJ mol\\S-1\\N nm\\S-1\\N)",
843 make_legend(outf,ngroups,isize0[0],index0[0],grpname,bCom,bMol,bDim,oenv);
847 outb = xvgropen(opt2fn("-ob",NFILE,fnm),"Box vector elements",
848 output_env_get_xvgr_tlabel(oenv),"(nm)",oenv);
850 xvgr_legend(outb,6,box_leg,oenv);
858 flags = flags | TRX_READ_V;
859 outt = xvgropen(opt2fn("-ot",NFILE,fnm),"Temperature",
860 output_env_get_xvgr_tlabel(oenv),"(K)", oenv);
861 make_legend(outt,ngroups,isize[0],index[0],grpname,bCom,bMol,bDum,oenv);
869 flags = flags | TRX_READ_V;
870 outekt = xvgropen(opt2fn("-ekt",NFILE,fnm),"Center of mass translation",
871 output_env_get_xvgr_tlabel(oenv),"Energy (kJ mol\\S-1\\N)",oenv);
872 make_legend(outekt,ngroups,isize[0],index[0],grpname,bCom,bMol,bDum,oenv);
880 flags = flags | TRX_READ_X | TRX_READ_V;
881 outekr = xvgropen(opt2fn("-ekr",NFILE,fnm),"Center of mass rotation",
882 output_env_get_xvgr_tlabel(oenv),"Energy (kJ mol\\S-1\\N)",oenv);
883 make_legend(outekr,ngroups,isize[0],index[0],grpname,bCom,bMol,bDum,oenv);
887 flags = flags | TRX_READ_V;
891 flags = flags | TRX_READ_X | TRX_READ_V;
895 flags = flags | TRX_READ_X | TRX_READ_F;
897 if ((flags == 0) && !bOB)
899 fprintf(stderr,"Please select one or more output file options\n");
903 read_first_frame(oenv,&status,ftp2fn(efTRX,NFILE,fnm),&fr,flags);
907 snew(sumxv,fr.natoms);
908 snew(sumv,fr.natoms);
912 snew(sumxf,fr.natoms);
913 snew(sumf,fr.natoms);
921 gpbc = gmx_rmpbc_init(&top.idef,ePBC,fr.natoms,fr.box);
926 time = output_env_conv_time(oenv,fr.time);
928 if (fr.bX && bNoJump && fr.bBox)
932 remove_jump(fr.box,fr.natoms,xp,fr.x);
938 for(i=0; i<fr.natoms; i++)
940 copy_rvec(fr.x[i],xp[i]);
944 if (fr.bX && bCom && bPBC)
946 gmx_rmpbc_trxfr(gpbc,&fr);
951 update_histo(isize[0],index[0],fr.v,&nvhisto,&vhisto,binwidth);
956 print_data(outx,time,fr.x,mass,bCom,ngroups,isize,index,bDim,sffmt);
963 frout.atoms = &top.atoms;
966 write_trx_x(status_out,&frout,mass,bCom,ngroups,isize,index);
970 print_data(outv,time,fr.v,mass,bCom,ngroups,isize,index,bDim,sffmt);
974 print_data(outf,time,fr.f,NULL,bCom,ngroups,isize,index,bDim,sffmt);
978 fprintf(outb,"\t%g",fr.time);
980 fr.box[XX][XX],fr.box[YY][YY],fr.box[ZZ][ZZ],
981 fr.box[YY][XX],fr.box[ZZ][XX],fr.box[ZZ][YY]);
986 fprintf(outt," %g",time);
987 for(i=0; i<ngroups; i++)
989 fprintf(outt,sffmt,temp(fr.v,mass,isize[i],index[i]));
995 fprintf(outekt," %g",time);
996 for(i=0; i<ngroups; i++)
998 fprintf(outekt,sffmt,ektrans(fr.v,mass,isize[i],index[i]));
1000 fprintf(outekt,"\n");
1002 if (bEKR && fr.bX && fr.bV)
1004 fprintf(outekr," %g",time);
1005 for(i=0; i<ngroups; i++)
1007 fprintf(outekr,sffmt,ekrot(fr.x,fr.v,mass,isize[i],index[i]));
1009 fprintf(outekr,"\n");
1015 for(i=0; i<fr.natoms; i++)
1017 rvec_inc(sumxv[i],fr.x[i]);
1023 for(i=0; i<fr.natoms; i++)
1024 rvec_inc(sumv[i],fr.v[i]);
1032 for(i=0; i<fr.natoms; i++)
1034 rvec_inc(sumxf[i],fr.x[i]);
1040 for(i=0; i<fr.natoms; i++)
1042 rvec_inc(sumf[i],fr.f[i]);
1048 } while(read_next_frame(oenv,status,&fr));
1050 gmx_rmpbc_done(gpbc);
1052 /* clean up a bit */
1055 if (bOX) ffclose(outx);
1056 if (bOXT) close_trx(status_out);
1057 if (bOV) ffclose(outv);
1058 if (bOF) ffclose(outf);
1059 if (bOB) ffclose(outb);
1060 if (bOT) ffclose(outt);
1061 if (bEKT) ffclose(outekt);
1062 if (bEKR) ffclose(outekr);
1066 print_histo(opt2fn("-vd",NFILE,fnm),nvhisto,vhisto,binwidth,oenv);
1069 if ((bCV || bCF) && (nr_vfr>1 || nr_ffr>1) && !bNoJump)
1071 fprintf(stderr,"WARNING: More than one frame was used for option -cv or -cf\n"
1072 "If atoms jump across the box you should use the -nojump option\n");
1077 write_pdb_bfac(opt2fn("-cv",NFILE,fnm),
1078 opt2fn("-av",NFILE,fnm),"average velocity",&(top.atoms),
1079 ePBC,topbox,isize[0],index[0],nr_xfr,sumxv,
1080 nr_vfr,sumv,bDim,scale,oenv);
1084 write_pdb_bfac(opt2fn("-cf",NFILE,fnm),
1085 opt2fn("-af",NFILE,fnm),"average force",&(top.atoms),
1086 ePBC,topbox,isize[0],index[0],nr_xfr,sumxf,
1087 nr_ffr,sumf,bDim,scale,oenv);
1091 view_all(oenv,NFILE, fnm);