3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Green Red Orange Magenta Azure Cyan Skyblue
62 #include "mtop_util.h"
76 real sumv, averv, maxv;
77 real *aver1, *aver2, *aver_3, *aver_6;
80 static void init5(int n)
86 static void reset5(void)
90 for (i = 0; (i < ntop); i++)
97 int tpcomp(const void *a, const void *b)
105 return (1e7*(tpb->v-tpa->v));
108 static void add5(int ndr, real viol)
113 for (i = 1; (i < ntop); i++)
115 if (top[i].v < top[mini].v)
120 if (viol > top[mini].v)
127 static void print5(FILE *fp)
131 qsort(top, ntop, sizeof(top[0]), tpcomp);
132 fprintf(fp, "Index:");
133 for (i = 0; (i < ntop); i++)
135 fprintf(fp, " %6d", top[i].n);
137 fprintf(fp, "\nViol: ");
138 for (i = 0; (i < ntop); i++)
140 fprintf(fp, " %6.3f", top[i].v);
145 static void check_viol(FILE *log, t_commrec *cr,
146 t_ilist *disres, t_iparams forceparams[],
148 t_forcerec *fr, t_pbc *pbc, t_graph *g, t_dr_result dr[],
149 int clust_id, int isize, atom_id index[], real vvindex[],
153 int i, j, nat, n, type, nviol, ndr, label;
154 real ener, rt, mviol, tviol, viol, lam, dvdl, drt;
155 static gmx_bool bFirst = TRUE;
167 forceatoms = disres->iatoms;
168 for (j = 0; (j < isize); j++)
172 nat = interaction_function[F_DISRES].nratoms+1;
173 for (i = 0; (i < disres->nr); )
175 type = forceatoms[i];
177 label = forceparams[type].disres.label;
180 fprintf(debug, "DISRE: ndr = %d, label = %d i=%d, n =%d\n",
185 gmx_fatal(FARGS, "tpr inconsistency. ndr = %d, label = %d\n", ndr, label);
191 while (((i+n) < disres->nr) &&
192 (forceparams[forceatoms[i+n]].disres.label == label));
194 calc_disres_R_6(cr->ms, n, &forceatoms[i], forceparams,
195 (const rvec*)x, pbc, fcd, NULL);
197 if (fcd->disres.Rt_6[0] <= 0)
199 gmx_fatal(FARGS, "ndr = %d, rt_6 = %f", ndr, fcd->disres.Rt_6[0]);
202 rt = pow(fcd->disres.Rt_6[0], -1.0/6.0);
203 dr[clust_id].aver1[ndr] += rt;
204 dr[clust_id].aver2[ndr] += sqr(rt);
206 dr[clust_id].aver_3[ndr] += drt;
207 dr[clust_id].aver_6[ndr] += fcd->disres.Rt_6[0];
209 ener = interaction_function[F_DISRES].ifunc(n, &forceatoms[i], forceparams,
210 (const rvec*)x, f, fr->fshift,
211 pbc, g, lam, &dvdl, NULL, fcd, NULL);
212 viol = fcd->disres.sumviol;
219 add5(forceparams[type].disres.label, viol);
226 for (j = 0; (j < isize); j++)
228 if (index[j] == forceparams[type].disres.label)
230 vvindex[j] = pow(fcd->disres.Rt_6[0], -1.0/6.0);
237 dr[clust_id].nv = nviol;
238 dr[clust_id].maxv = mviol;
239 dr[clust_id].sumv = tviol;
240 dr[clust_id].averv = tviol/ndr;
241 dr[clust_id].nframes++;
245 fprintf(stderr, "\nThere are %d restraints and %d pairs\n",
246 ndr, disres->nr/nat);
258 real up1, r, rT3, rT6, viol, violT3, violT6;
261 static int drs_comp(const void *a, const void *b)
265 da = (t_dr_stats *)a;
266 db = (t_dr_stats *)b;
268 if (da->viol > db->viol)
272 else if (da->viol < db->viol)
282 static void dump_dump(FILE *log, int ndr, t_dr_stats drs[])
284 static const char *core[] = { "All restraints", "Core restraints" };
285 static const char *tp[] = { "linear", "third power", "sixth power" };
286 real viol_tot, viol_max, viol = 0;
291 for (bCore = FALSE; (bCore <= TRUE); bCore++)
293 for (kkk = 0; (kkk < 3); kkk++)
299 for (i = 0; (i < ndr); i++)
301 if (!bCore || (bCore && drs[i].bCore))
309 viol = drs[i].violT3;
312 viol = drs[i].violT6;
315 gmx_incons("Dumping violations");
317 viol_max = max(viol_max, viol);
326 if ((nrestr > 0) || (bCore && (nrestr < ndr)))
329 fprintf(log, "+++++++ %s ++++++++\n", core[bCore]);
330 fprintf(log, "+++++++ Using %s averaging: ++++++++\n", tp[kkk]);
331 fprintf(log, "Sum of violations: %8.3f nm\n", viol_tot);
334 fprintf(log, "Average violation: %8.3f nm\n", viol_tot/nrestr);
336 fprintf(log, "Largest violation: %8.3f nm\n", viol_max);
337 fprintf(log, "Number of violated restraints: %d/%d\n", nviol, nrestr);
343 static void dump_viol(FILE *log, int ndr, t_dr_stats *drs, gmx_bool bLinear)
347 fprintf(log, "Restr. Core Up1 <r> <rT3> <rT6> <viol><violT3><violT6>\n");
348 for (i = 0; (i < ndr); i++)
350 if (bLinear && (drs[i].viol == 0))
354 fprintf(log, "%6d%5s%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f\n",
355 drs[i].index, yesno_names[drs[i].bCore],
356 drs[i].up1, drs[i].r, drs[i].rT3, drs[i].rT6,
357 drs[i].viol, drs[i].violT3, drs[i].violT6);
361 static gmx_bool is_core(int i, int isize, atom_id index[])
364 gmx_bool bIC = FALSE;
366 for (kk = 0; !bIC && (kk < isize); kk++)
368 bIC = (index[kk] == i);
374 static void dump_stats(FILE *log, int nsteps, int ndr, t_ilist *disres,
375 t_iparams ip[], t_dr_result *dr,
376 int isize, atom_id index[], t_atoms *atoms)
382 fprintf(log, "++++++++++++++ STATISTICS ++++++++++++++++++++++++\n");
385 nra = interaction_function[F_DISRES].nratoms+1;
386 for (i = 0; (i < ndr); i++)
388 /* Search for the appropriate restraint */
389 for (; (j < disres->nr); j += nra)
391 if (ip[disres->iatoms[j]].disres.label == i)
397 drs[i].bCore = is_core(i, isize, index);
398 drs[i].up1 = ip[disres->iatoms[j]].disres.up1;
399 drs[i].r = dr->aver1[i]/nsteps;
400 drs[i].rT3 = pow(dr->aver_3[i]/nsteps, -1.0/3.0);
401 drs[i].rT6 = pow(dr->aver_6[i]/nsteps, -1.0/6.0);
402 drs[i].viol = max(0, drs[i].r-drs[i].up1);
403 drs[i].violT3 = max(0, drs[i].rT3-drs[i].up1);
404 drs[i].violT6 = max(0, drs[i].rT6-drs[i].up1);
407 int j1 = disres->iatoms[j+1];
408 int j2 = disres->iatoms[j+2];
409 atoms->pdbinfo[j1].bfac += drs[i].violT3*5;
410 atoms->pdbinfo[j2].bfac += drs[i].violT3*5;
413 dump_viol(log, ndr, drs, FALSE);
415 fprintf(log, "+++ Sorted by linear averaged violations: +++\n");
416 qsort(drs, ndr, sizeof(drs[0]), drs_comp);
417 dump_viol(log, ndr, drs, TRUE);
419 dump_dump(log, ndr, drs);
424 static void dump_clust_stats(FILE *fp, int ndr, t_ilist *disres,
425 t_iparams ip[], t_blocka *clust, t_dr_result dr[],
426 char *clust_name[], int isize, atom_id index[])
428 int i, j, k, nra, mmm = 0;
429 double sumV, maxV, sumVT3, sumVT6, maxVT3, maxVT6;
433 fprintf(fp, "++++++++++++++ STATISTICS ++++++++++++++++++++++\n");
434 fprintf(fp, "Cluster NFrames SumV MaxV SumVT MaxVT SumVS MaxVS\n");
438 for (k = 0; (k < clust->nr); k++)
440 if (dr[k].nframes == 0)
444 if (dr[k].nframes != (clust->index[k+1]-clust->index[k]))
446 gmx_fatal(FARGS, "Inconsistency in cluster %s.\n"
447 "Found %d frames in trajectory rather than the expected %d\n",
448 clust_name[k], dr[k].nframes,
449 clust->index[k+1]-clust->index[k]);
453 gmx_fatal(FARGS, "Inconsistency with cluster %d. Invalid name", k);
456 nra = interaction_function[F_DISRES].nratoms+1;
457 sumV = sumVT3 = sumVT6 = maxV = maxVT3 = maxVT6 = 0;
458 for (i = 0; (i < ndr); i++)
460 /* Search for the appropriate restraint */
461 for (; (j < disres->nr); j += nra)
463 if (ip[disres->iatoms[j]].disres.label == i)
469 drs[i].bCore = is_core(i, isize, index);
470 drs[i].up1 = ip[disres->iatoms[j]].disres.up1;
471 drs[i].r = dr[k].aver1[i]/dr[k].nframes;
472 if ((dr[k].aver_3[i] <= 0) || (dr[k].aver_3[i] != dr[k].aver_3[i]))
474 gmx_fatal(FARGS, "dr[%d].aver_3[%d] = %f", k, i, dr[k].aver_3[i]);
476 drs[i].rT3 = pow(dr[k].aver_3[i]/dr[k].nframes, -1.0/3.0);
477 drs[i].rT6 = pow(dr[k].aver_6[i]/dr[k].nframes, -1.0/6.0);
478 drs[i].viol = max(0, drs[i].r-drs[i].up1);
479 drs[i].violT3 = max(0, drs[i].rT3-drs[i].up1);
480 drs[i].violT6 = max(0, drs[i].rT6-drs[i].up1);
482 sumVT3 += drs[i].violT3;
483 sumVT6 += drs[i].violT6;
484 maxV = max(maxV, drs[i].viol);
485 maxVT3 = max(maxVT3, drs[i].violT3);
486 maxVT6 = max(maxVT6, drs[i].violT6);
488 if (strcmp(clust_name[k], "1000") == 0)
492 fprintf(fp, "%-10s%6d%8.3f %8.3f %8.3f %8.3f %8.3f %8.3f\n",
494 dr[k].nframes, sumV, maxV, sumVT3, maxVT3, sumVT6, maxVT6);
501 static void init_dr_res(t_dr_result *dr, int ndr)
503 snew(dr->aver1, ndr+1);
504 snew(dr->aver2, ndr+1);
505 snew(dr->aver_3, ndr+1);
506 snew(dr->aver_6, ndr+1);
514 static void dump_disre_matrix(const char *fn, t_dr_result *dr, int ndr,
515 int nsteps, t_idef *idef, gmx_mtop_t *mtop,
516 real max_dr, int nlevels, gmx_bool bThird)
520 int n_res, a_offset, mb, mol, a;
522 int iii, i, j, nra, nratoms, tp, ri, rj, index, nlabel, label;
523 atom_id ai, aj, *ptr;
524 real **matrix, *t_res, hi, *w_dr, rav, rviol;
525 t_rgb rlo = { 1, 1, 1 };
526 t_rgb rhi = { 0, 0, 0 };
532 snew(resnr, mtop->natoms);
535 for (mb = 0; mb < mtop->nmolblock; mb++)
537 atoms = &mtop->moltype[mtop->molblock[mb].type].atoms;
538 for (mol = 0; mol < mtop->molblock[mb].nmol; mol++)
540 for (a = 0; a < atoms->nr; a++)
542 resnr[a_offset+a] = n_res + atoms->atom[a].resind;
544 n_res += atoms->nres;
545 a_offset += atoms->nr;
550 for (i = 0; (i < n_res); i++)
555 for (i = 0; (i < n_res); i++)
557 snew(matrix[i], n_res);
559 nratoms = interaction_function[F_DISRES].nratoms;
560 nra = (idef->il[F_DISRES].nr/(nratoms+1));
566 for (i = 0; (i < idef->il[F_DISRES].nr); i += nratoms+1)
568 tp = idef->il[F_DISRES].iatoms[i];
569 label = idef->iparams[tp].disres.label;
573 /* Set index pointer */
577 gmx_fatal(FARGS, "nlabel is %d, label = %d", nlabel, label);
581 gmx_fatal(FARGS, "ndr = %d, index = %d", ndr, index);
583 /* Update the weight */
584 w_dr[index] = 1.0/nlabel;
593 printf("nlabel = %d, index = %d, ndr = %d\n", nlabel, index, ndr);
595 for (i = 0; (i < ndr); i++)
597 for (j = ptr[i]; (j < ptr[i+1]); j += nratoms+1)
599 tp = idef->il[F_DISRES].iatoms[j];
600 ai = idef->il[F_DISRES].iatoms[j+1];
601 aj = idef->il[F_DISRES].iatoms[j+2];
607 rav = pow(dr->aver_3[i]/nsteps, -1.0/3.0);
611 rav = dr->aver1[i]/nsteps;
615 fprintf(debug, "DR %d, atoms %d, %d, distance %g\n", i, ai, aj, rav);
617 rviol = max(0, rav-idef->iparams[tp].disres.up1);
618 matrix[ri][rj] += w_dr[i]*rviol;
619 matrix[rj][ri] += w_dr[i]*rviol;
620 hi = max(hi, matrix[ri][rj]);
621 hi = max(hi, matrix[rj][ri]);
631 printf("Warning: the maxdr that you have specified (%g) is smaller than\nthe largest value in your simulation (%g)\n", max_dr, hi);
635 printf("Highest level in the matrix will be %g\n", hi);
636 fp = ffopen(fn, "w");
637 write_xpm(fp, 0, "Distance Violations", "<V> (nm)", "Residue", "Residue",
638 n_res, n_res, t_res, t_res, matrix, 0, hi, rlo, rhi, &nlevels);
642 int gmx_disre(int argc, char *argv[])
644 const char *desc[] = {
645 "[TT]g_disre[tt] computes violations of distance restraints.",
646 "If necessary, all protons can be added to a protein molecule ",
647 "using the [TT]g_protonate[tt] program.[PAR]",
648 "The program always",
649 "computes the instantaneous violations rather than time-averaged,",
650 "because this analysis is done from a trajectory file afterwards",
651 "it does not make sense to use time averaging. However,",
652 "the time averaged values per restraint are given in the log file.[PAR]",
653 "An index file may be used to select specific restraints for",
655 "When the optional [TT]-q[tt] flag is given a [TT].pdb[tt] file coloured by the",
656 "amount of average violations.[PAR]",
657 "When the [TT]-c[tt] option is given, an index file will be read",
658 "containing the frames in your trajectory corresponding to the clusters",
659 "(defined in another manner) that you want to analyze. For these clusters",
660 "the program will compute average violations using the third power",
661 "averaging algorithm and print them in the log file."
664 static int nlevels = 20;
665 static real max_dr = 0;
666 static gmx_bool bThird = TRUE;
668 { "-ntop", FALSE, etINT, {&ntop},
669 "Number of large violations that are stored in the log file every step" },
670 { "-maxdr", FALSE, etREAL, {&max_dr},
671 "Maximum distance violation in matrix output. If less than or equal to 0 the maximum will be determined by the data." },
672 { "-nlevels", FALSE, etINT, {&nlevels},
673 "Number of levels in the matrix output" },
674 { "-third", FALSE, etBOOL, {&bThird},
675 "Use inverse third power averaging or linear for matrix output" }
678 FILE *out = NULL, *aver = NULL, *numv = NULL, *maxxv = NULL, *xvg = NULL;
684 t_atoms *atoms = NULL;
690 int ntopatoms, natoms, i, j, kkk;
693 rvec *x, *f, *xav = NULL;
697 atom_id *index = NULL, *ind_fit = NULL;
699 t_cluster_ndx *clust = NULL;
700 t_dr_result dr, *dr_clust = NULL;
702 real *vvindex = NULL, *w_rls = NULL;
704 t_pbc pbc, *pbc_null;
708 gmx_rmpbc_t gpbc = NULL;
711 { efTPX, NULL, NULL, ffREAD },
712 { efTRX, "-f", NULL, ffREAD },
713 { efXVG, "-ds", "drsum", ffWRITE },
714 { efXVG, "-da", "draver", ffWRITE },
715 { efXVG, "-dn", "drnum", ffWRITE },
716 { efXVG, "-dm", "drmax", ffWRITE },
717 { efXVG, "-dr", "restr", ffWRITE },
718 { efLOG, "-l", "disres", ffWRITE },
719 { efNDX, NULL, "viol", ffOPTRD },
720 { efPDB, "-q", "viol", ffOPTWR },
721 { efNDX, "-c", "clust", ffOPTRD },
722 { efXPM, "-x", "matrix", ffOPTWR }
724 #define NFILE asize(fnm)
726 cr = init_par(&argc, &argv);
727 parse_common_args(&argc, argv, PCA_CAN_TIME | PCA_CAN_VIEW | PCA_BE_NICE,
728 NFILE, fnm, asize(pa), pa, asize(desc), desc, 0, NULL, &oenv);
730 gmx_log_open(ftp2fn(efLOG, NFILE, fnm), cr, FALSE, 0, &fplog);
737 read_tpxheader(ftp2fn(efTPX, NFILE, fnm), &header, FALSE, NULL, NULL);
738 snew(xtop, header.natoms);
739 read_tpx(ftp2fn(efTPX, NFILE, fnm), &ir, box, &ntopatoms, xtop, NULL, NULL, &mtop);
740 bPDB = opt2bSet("-q", NFILE, fnm);
743 snew(xav, ntopatoms);
744 snew(ind_fit, ntopatoms);
745 snew(w_rls, ntopatoms);
746 for (kkk = 0; (kkk < ntopatoms); kkk++)
753 *atoms = gmx_mtop_global_atoms(&mtop);
755 if (atoms->pdbinfo == NULL)
757 snew(atoms->pdbinfo, atoms->nr);
761 top = gmx_mtop_generate_local_top(&mtop, &ir);
765 if (ir.ePBC != epbcNONE)
767 if (ir.bPeriodicMols)
773 g = mk_graph(fplog, &top->idef, 0, mtop.natoms, FALSE, FALSE);
777 if (ftp2bSet(efNDX, NFILE, fnm))
779 rd_index(ftp2fn(efNDX, NFILE, fnm), 1, &isize, &index, &grpname);
780 xvg = xvgropen(opt2fn("-dr", NFILE, fnm), "Individual Restraints", "Time (ps)",
782 snew(vvindex, isize);
784 for (i = 0; (i < isize); i++)
788 sprintf(leg[i], "index %d", index[i]);
790 xvgr_legend(xvg, isize, (const char**)leg, oenv);
798 init_disres(fplog, &mtop, &ir, NULL, FALSE, &fcd, NULL, FALSE);
800 natoms = read_first_x(oenv, &status, ftp2fn(efTRX, NFILE, fnm), &t, &x, box);
803 init_dr_res(&dr, fcd.disres.nres);
804 if (opt2bSet("-c", NFILE, fnm))
806 clust = cluster_index(fplog, opt2fn("-c", NFILE, fnm));
807 snew(dr_clust, clust->clust->nr+1);
808 for (i = 0; (i <= clust->clust->nr); i++)
810 init_dr_res(&dr_clust[i], fcd.disres.nres);
815 out = xvgropen(opt2fn("-ds", NFILE, fnm),
816 "Sum of Violations", "Time (ps)", "nm", oenv);
817 aver = xvgropen(opt2fn("-da", NFILE, fnm),
818 "Average Violation", "Time (ps)", "nm", oenv);
819 numv = xvgropen(opt2fn("-dn", NFILE, fnm),
820 "# Violations", "Time (ps)", "#", oenv);
821 maxxv = xvgropen(opt2fn("-dm", NFILE, fnm),
822 "Largest Violation", "Time (ps)", "nm", oenv);
825 mdatoms = init_mdatoms(fplog, &mtop, ir.efep != efepNO);
826 atoms2md(&mtop, &ir, 0, NULL, 0, mtop.natoms, mdatoms);
827 update_mdatoms(mdatoms, ir.fepvals->init_lambda);
829 fprintf(fplog, "Made forcerec\n");
830 init_forcerec(fplog, oenv, fr, NULL, &ir, &mtop, cr, box, FALSE,
831 NULL, NULL, NULL, NULL, NULL, FALSE, -1);
833 if (ir.ePBC != epbcNONE)
835 gpbc = gmx_rmpbc_init(&top->idef, ir.ePBC, natoms, box);
841 if (ir.ePBC != epbcNONE)
843 if (ir.bPeriodicMols)
845 set_pbc(&pbc, ir.ePBC, box);
849 gmx_rmpbc(gpbc, natoms, box, x);
855 if (j > clust->maxframe)
857 gmx_fatal(FARGS, "There are more frames in the trajectory than in the cluster index file. t = %8f\n", t);
859 my_clust = clust->inv_clust[j];
860 range_check(my_clust, 0, clust->clust->nr);
861 check_viol(fplog, cr, &(top->idef.il[F_DISRES]),
863 x, f, fr, pbc_null, g, dr_clust, my_clust, isize, index, vvindex, &fcd);
867 check_viol(fplog, cr, &(top->idef.il[F_DISRES]),
869 x, f, fr, pbc_null, g, &dr, 0, isize, index, vvindex, &fcd);
873 reset_x(atoms->nr, ind_fit, atoms->nr, NULL, x, w_rls);
874 do_fit(atoms->nr, w_rls, x, x);
877 /* Store the first frame of the trajectory as 'characteristic'
878 * for colouring with violations.
880 for (kkk = 0; (kkk < atoms->nr); kkk++)
882 copy_rvec(x[kkk], xav[kkk]);
890 fprintf(xvg, "%10g", t);
891 for (i = 0; (i < isize); i++)
893 fprintf(xvg, " %10g", vvindex[i]);
897 fprintf(out, "%10g %10g\n", t, dr.sumv);
898 fprintf(aver, "%10g %10g\n", t, dr.averv);
899 fprintf(maxxv, "%10g %10g\n", t, dr.maxv);
900 fprintf(numv, "%10g %10d\n", t, dr.nv);
904 while (read_next_x(oenv, status, &t, natoms, x, box));
906 if (ir.ePBC != epbcNONE)
908 gmx_rmpbc_done(gpbc);
913 dump_clust_stats(fplog, fcd.disres.nres, &(top->idef.il[F_DISRES]),
914 top->idef.iparams, clust->clust, dr_clust,
915 clust->grpname, isize, index);
919 dump_stats(fplog, j, fcd.disres.nres, &(top->idef.il[F_DISRES]),
920 top->idef.iparams, &dr, isize, index,
921 bPDB ? atoms : NULL);
924 write_sto_conf(opt2fn("-q", NFILE, fnm),
925 "Coloured by average violation in Angstrom",
926 atoms, xav, NULL, ir.ePBC, box);
928 dump_disre_matrix(opt2fn_null("-x", NFILE, fnm), &dr, fcd.disres.nres,
929 j, &top->idef, &mtop, max_dr, nlevels, bThird);
937 do_view(oenv, opt2fn("-dr", NFILE, fnm), "-nxy");
939 do_view(oenv, opt2fn("-dn", NFILE, fnm), "-nxy");
940 do_view(oenv, opt2fn("-da", NFILE, fnm), "-nxy");
941 do_view(oenv, opt2fn("-ds", NFILE, fnm), "-nxy");
942 do_view(oenv, opt2fn("-dm", NFILE, fnm), "-nxy");
947 gmx_log_close(fplog);