-static void do_ballistic(const char *balFile, int nData,
- real *t, real **val, int nSet,
- real balTime, int nBalExp,
- const output_env_t oenv)
-{
- double **ctd = NULL, *td = NULL;
- t_gemParams *GP = init_gemParams(0, 0, t, nData, 0, 0, 0, balTime, nBalExp);
- static char *leg[] = {"Ac'(t)"};
- FILE *fp;
- int i, set;
-
- if (GP->ballistic/GP->tDelta >= GP->nExpFit*2+1)
- {
- snew(ctd, nSet);
- snew(td, nData);
-
- fp = xvgropen(balFile, "Hydrogen Bond Autocorrelation", "Time (ps)", "C'(t)", oenv);
- xvgr_legend(fp, asize(leg), (const char**)leg, oenv);
-
- for (set = 0; set < nSet; set++)
- {
- snew(ctd[set], nData);
- for (i = 0; i < nData; i++)
- {
- ctd[set][i] = (double)val[set][i];
- if (set == 0)
- {
- td[i] = (double)t[i];
- }
- }
-
- takeAwayBallistic(ctd[set], td, nData, GP->ballistic, GP->nExpFit, GP->bDt);
- }
-
- for (i = 0; i < nData; i++)
- {
- fprintf(fp, " %g", t[i]);
- for (set = 0; set < nSet; set++)
- {
- fprintf(fp, " %g", ctd[set][i]);
- }
- fprintf(fp, "\n");
- }
-
-
- for (set = 0; set < nSet; set++)
- {
- sfree(ctd[set]);
- }
- sfree(ctd);
- sfree(td);
- xvgrclose(fp);
- }
- else
- {
- printf("Number of data points is less than the number of parameters to fit\n."
- "The system is underdetermined, hence no ballistic term can be found.\n\n");
- }
-}
-
-static void do_geminate(const char *gemFile, int nData,
- real *t, real **val, int nSet,
- const real D, const real rcut, const real balTime,
- const int nFitPoints, const real begFit, const real endFit,
- const output_env_t oenv)
-{
- double **ctd = NULL, **ctdGem = NULL, *td = NULL;
- t_gemParams *GP = init_gemParams(rcut, D, t, nData, nFitPoints,
- begFit, endFit, balTime, 1);
- const char *leg[] = {"Ac\\sgem\\N(t)"};
- FILE *fp;
- int i, set;
-
- snew(ctd, nSet);
- snew(ctdGem, nSet);
- snew(td, nData);
-
- fp = xvgropen(gemFile, "Hydrogen Bond Autocorrelation", "Time (ps)", "C'(t)", oenv);
- xvgr_legend(fp, asize(leg), leg, oenv);
-
- for (set = 0; set < nSet; set++)
- {
- snew(ctd[set], nData);
- snew(ctdGem[set], nData);
- for (i = 0; i < nData; i++)
- {
- ctd[set][i] = (double)val[set][i];
- if (set == 0)
- {
- td[i] = (double)t[i];
- }
- }
- fitGemRecomb(ctd[set], td, &(ctd[set]), nData, GP);
- }
-
- for (i = 0; i < nData; i++)
- {
- fprintf(fp, " %g", t[i]);
- for (set = 0; set < nSet; set++)
- {
- fprintf(fp, " %g", ctdGem[set][i]);
- }
- fprintf(fp, "\n");
- }
-
- for (set = 0; set < nSet; set++)
- {
- sfree(ctd[set]);
- sfree(ctdGem[set]);
- }
- sfree(ctd);
- sfree(ctdGem);
- sfree(td);
- xvgrclose(fp);
-}
-