src/gromacs/analysisdata/modules/plot.cpp

   1 /*
   2  *
   3  *                This source code is part of
   4  *
   5  *                 G   R   O   M   A   C   S
   6  *
   7  *          GROningen MAchine for Chemical Simulations
   8  *
   9  * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
  10  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  11  * Copyright (c) 2001-2009, The GROMACS development team,
  12  * check out http://www.gromacs.org for more information.
  13
  14  * This program is free software; you can redistribute it and/or
  15  * modify it under the terms of the GNU General Public License
  16  * as published by the Free Software Foundation; either version 2
  17  * of the License, or (at your option) any later version.
  18  *
  19  * If you want to redistribute modifications, please consider that
  20  * scientific software is very special. Version control is crucial -
  21  * bugs must be traceable. We will be happy to consider code for
  22  * inclusion in the official distribution, but derived work must not
  23  * be called official GROMACS. Details are found in the README & COPYING
  24  * files - if they are missing, get the official version at www.gromacs.org.
  25  *
  26  * To help us fund GROMACS development, we humbly ask that you cite
  27  * the papers on the package - you can find them in the top README file.
  28  *
  29  * For more info, check our website at http://www.gromacs.org
  30  */
  31 /*! \internal \file
  32  * \brief
  33  * Implements classes in plot.h.
  34  *
  35  * \ingroup module_analysisdata
  36  * \author Teemu Murtola <teemu.murtola@cbr.su.se>
  37  */
  38 #include "gromacs/analysisdata/modules/plot.h"
  39
  40 #include <string>
  41 #include <vector>
  42
  43 #include <cstdio>
  44 #include <cstring>
  45
  46 #include <boost/shared_ptr.hpp>
  47
  48 #include "gromacs/legacyheaders/gmxfio.h"
  49 #include "gromacs/legacyheaders/oenv.h"
  50 #include "gromacs/legacyheaders/vec.h"
  51 #include "gromacs/legacyheaders/xvgr.h"
  52
  53 #include "gromacs/options/basicoptions.h"
  54 #include "gromacs/analysisdata/dataframe.h"
  55 #include "gromacs/options/options.h"
  56 #include "gromacs/options/timeunitmanager.h"
  57 #include "gromacs/selection/selectioncollection.h"
  58 #include "gromacs/utility/exceptions.h"
  59 #include "gromacs/utility/gmxassert.h"
  60 #include "gromacs/utility/stringutil.h"
  61
  62 namespace
  63 {
  64
  65 //! Enum values for plot formats.
  66 const char *const g_plotFormats[] = {
  67     "none", "xmgrace", "xmgr", NULL
  68 };
  69
  70 } // namespace
  71
  72 namespace gmx
  73 {
  74
  75 /********************************************************************
  76  * AnalysisDataPlotSettings
  77  */
  78
  79 AnalysisDataPlotSettings::AnalysisDataPlotSettings()
  80     : selections_(NULL), timeUnit_(eTimeUnit_ps), plotFormat_(1)
  81 {
  82 }
  83
  84 void
  85 AnalysisDataPlotSettings::setSelectionCollection(const SelectionCollection *selections)
  86 {
  87     selections_ = selections;
  88 }
  89
  90
  91 void
  92 AnalysisDataPlotSettings::addOptions(Options *options)
  93 {
  94     options->addOption(StringOption("xvg").enumValue(g_plotFormats)
  95                            .defaultValue("xmgrace")
  96                            .storeEnumIndex(&plotFormat_)
  97                            .description("Plot formatting"));
  98 }
  99
 100
 101 /********************************************************************
 102  * AbstractPlotModule::Impl
 103  */
 104
 105 class AbstractPlotModule::Impl
 106 {
 107     public:
 108         explicit Impl(const AnalysisDataPlotSettings &settings);
 109         ~Impl();
 110
 111         void closeFile();
 112
 113         AnalysisDataPlotSettings  settings_;
 114         std::string               filename_;
 115         FILE                     *fp_;
 116
 117         bool                      bPlain_;
 118         bool                      gOmitX_;
 119         std::string               title_;
 120         std::string               subtitle_;
 121         std::string               xlabel_;
 122         std::string               ylabel_;
 123         std::vector<std::string>  legend_;
 124         char                      xformat_[15];
 125         char                      yformat_[15];
 126         real                      xscale_;
 127 };
 128
 129 AbstractPlotModule::Impl::Impl(const AnalysisDataPlotSettings &settings)
 130     : settings_(settings), fp_(NULL), bPlain_(false), gOmitX_(false),
 131       xscale_(1.0)
 132 {
 133     strcpy(xformat_, "%11.3f");
 134     strcpy(yformat_, " %8.3f");
 135 }
 136
 137 AbstractPlotModule::Impl::~Impl()
 138 {
 139     closeFile();
 140 }
 141
 142
 143 void
 144 AbstractPlotModule::Impl::closeFile()
 145 {
 146     if (fp_ != NULL)
 147     {
 148         if (bPlain_)
 149         {
 150             gmx_fio_fclose(fp_);
 151         }
 152         else
 153         {
 154             xvgrclose(fp_);
 155         }
 156         fp_ = NULL;
 157     }
 158 }
 159
 160
 161 /********************************************************************
 162  * AbstractPlotModule
 163  */
 164 /*! \cond libapi */
 165 AbstractPlotModule::AbstractPlotModule()
 166     : impl_(new Impl(AnalysisDataPlotSettings()))
 167 {
 168 }
 169
 170 AbstractPlotModule::AbstractPlotModule(const AnalysisDataPlotSettings &settings)
 171     : impl_(new Impl(settings))
 172 {
 173 }
 174 //! \endcond
 175
 176 AbstractPlotModule::~AbstractPlotModule()
 177 {
 178 }
 179
 180
 181 void
 182 AbstractPlotModule::setSettings(const AnalysisDataPlotSettings &settings)
 183 {
 184     impl_->settings_ = settings;
 185 }
 186
 187
 188 void
 189 AbstractPlotModule::setFileName(const std::string &filename)
 190 {
 191     impl_->filename_ = filename;
 192 }
 193
 194
 195 void
 196 AbstractPlotModule::setPlainOutput(bool bPlain)
 197 {
 198     impl_->bPlain_ = bPlain;
 199 }
 200
 201
 202 void
 203 AbstractPlotModule::setOmitX(bool bOmitX)
 204 {
 205     impl_->gOmitX_ = bOmitX;
 206 }
 207
 208
 209 void
 210 AbstractPlotModule::setTitle(const char *title)
 211 {
 212     impl_->title_ = title;
 213 }
 214
 215
 216 void
 217 AbstractPlotModule::setSubtitle(const char *subtitle)
 218 {
 219     impl_->subtitle_ = subtitle;
 220 }
 221
 222
 223 void
 224 AbstractPlotModule::setXLabel(const char *label)
 225 {
 226     impl_->xlabel_ = label;
 227 }
 228
 229
 230 void
 231 AbstractPlotModule::setXAxisIsTime()
 232 {
 233     TimeUnitManager manager(impl_->settings_.timeUnit());
 234     impl_->xlabel_ = formatString("Time (%s)", manager.timeUnitAsString());
 235     impl_->xscale_ = manager.inverseTimeScaleFactor();
 236 }
 237
 238
 239 void
 240 AbstractPlotModule::setYLabel(const char *label)
 241 {
 242     impl_->ylabel_ = label;
 243 }
 244
 245
 246 void
 247 AbstractPlotModule::setLegend(int nsets, const char * const *setname)
 248 {
 249     impl_->legend_.reserve(impl_->legend_.size() + nsets);
 250     for (int i = 0; i < nsets; ++i)
 251     {
 252         appendLegend(setname[i]);
 253     }
 254 }
 255
 256
 257 void
 258 AbstractPlotModule::appendLegend(const char *setname)
 259 {
 260     impl_->legend_.push_back(setname);
 261 }
 262
 263
 264 void
 265 AbstractPlotModule::setXFormat(int width, int precision, char format)
 266 {
 267     GMX_RELEASE_ASSERT(width >= 0 && precision >= 0
 268                        && width <= 99 && precision <= 99,
 269                        "Invalid width or precision");
 270     GMX_RELEASE_ASSERT(strchr("eEfFgG", format) != NULL,
 271                        "Invalid format specifier");
 272     sprintf(impl_->xformat_, "%%%d.%d%c", width, precision, format);
 273 }
 274
 275
 276 void
 277 AbstractPlotModule::setYFormat(int width, int precision, char format)
 278 {
 279     GMX_RELEASE_ASSERT(width >= 0 && precision >= 0
 280                        && width <= 99 && precision <= 99,
 281                        "Invalid width or precision");
 282     GMX_RELEASE_ASSERT(strchr("eEfFgG", format) != NULL,
 283                        "Invalid format specifier");
 284     sprintf(impl_->yformat_, " %%%d.%d%c", width, precision, format);
 285 }
 286
 287
 288 int
 289 AbstractPlotModule::flags() const
 290 {
 291     return efAllowMulticolumn | efAllowMultipoint;
 292 }
 293
 294
 295 void
 296 AbstractPlotModule::dataStarted(AbstractAnalysisData *data)
 297 {
 298     if (!impl_->filename_.empty())
 299     {
 300         if (impl_->bPlain_)
 301         {
 302             impl_->fp_ = gmx_fio_fopen(impl_->filename_.c_str(), "w");
 303         }
 304         else
 305         {
 306             time_unit_t  time_unit
 307                 = static_cast<time_unit_t>(impl_->settings_.timeUnit() + 1);
 308             xvg_format_t xvg_format
 309                 = (impl_->settings_.plotFormat() > 0
 310                    ? static_cast<xvg_format_t>(impl_->settings_.plotFormat())
 311                    : exvgNONE);
 312             output_env_t                  oenv;
 313             output_env_init(&oenv, 0, NULL, time_unit, FALSE, xvg_format, 0, 0);
 314             boost::shared_ptr<output_env> oenvGuard(oenv, &output_env_done);
 315             impl_->fp_ = xvgropen(impl_->filename_.c_str(), impl_->title_.c_str(),
 316                                   impl_->xlabel_.c_str(), impl_->ylabel_.c_str(),
 317                                   oenv);
 318             const SelectionCollection *selections
 319                 = impl_->settings_.selectionCollection();
 320             if (selections != NULL)
 321             {
 322                 selections->printXvgrInfo(impl_->fp_, oenv);
 323             }
 324             if (!impl_->subtitle_.empty())
 325             {
 326                 xvgr_subtitle(impl_->fp_, impl_->subtitle_.c_str(), oenv);
 327             }
 328             if (output_env_get_print_xvgr_codes(oenv)
 329                 && !impl_->legend_.empty())
 330             {
 331                 std::vector<const char *> legend;
 332                 legend.reserve(impl_->legend_.size());
 333                 for (size_t i = 0; i < impl_->legend_.size(); ++i)
 334                 {
 335                     legend.push_back(impl_->legend_[i].c_str());
 336                 }
 337                 xvgr_legend(impl_->fp_, legend.size(), &legend[0], oenv);
 338             }
 339         }
 340     }
 341 }
 342
 343
 344 void
 345 AbstractPlotModule::frameStarted(const AnalysisDataFrameHeader &frame)
 346 {
 347     if (!isFileOpen())
 348     {
 349         return;
 350     }
 351     if (!impl_->gOmitX_)
 352     {
 353         std::fprintf(impl_->fp_, impl_->xformat_, frame.x() * impl_->xscale_);
 354     }
 355 }
 356
 357
 358 void
 359 AbstractPlotModule::frameFinished(const AnalysisDataFrameHeader & /*header*/)
 360 {
 361     if (!isFileOpen())
 362     {
 363         return;
 364     }
 365     std::fprintf(impl_->fp_, "\n");
 366 }
 367
 368
 369 void
 370 AbstractPlotModule::dataFinished()
 371 {
 372     impl_->closeFile();
 373 }
 374
 375 /*! \cond libapi */
 376 bool
 377 AbstractPlotModule::isFileOpen() const
 378 {
 379     return impl_->fp_ != NULL;
 380 }
 381
 382
 383 void
 384 AbstractPlotModule::writeValue(real value) const
 385 {
 386     GMX_ASSERT(isFileOpen(), "File not opened, but write attempted");
 387     std::fprintf(impl_->fp_, impl_->yformat_, value);
 388 }
 389 //! \endcond
 390
 391 /********************************************************************
 392  * DataPlotModule
 393  */
 394
 395 AnalysisDataPlotModule::AnalysisDataPlotModule()
 396 {
 397 }
 398
 399 AnalysisDataPlotModule::AnalysisDataPlotModule(
 400         const AnalysisDataPlotSettings &settings)
 401     : AbstractPlotModule(settings)
 402 {
 403 }
 404
 405
 406 void
 407 AnalysisDataPlotModule::pointsAdded(const AnalysisDataPointSetRef &points)
 408 {
 409     if (!isFileOpen())
 410     {
 411         return;
 412     }
 413     for (int i = 0; i < points.columnCount(); ++i)
 414     {
 415         writeValue(points.y(i));
 416     }
 417 }
 418
 419
 420 /********************************************************************
 421  * DataVectorPlotModule
 422  */
 423
 424 AnalysisDataVectorPlotModule::AnalysisDataVectorPlotModule()
 425 {
 426     for (int i = 0; i < DIM; ++i)
 427     {
 428         bWrite_[i] = true;
 429     }
 430     bWrite_[DIM] = false;
 431 }
 432
 433
 434 AnalysisDataVectorPlotModule::AnalysisDataVectorPlotModule(
 435         const AnalysisDataPlotSettings &settings)
 436     : AbstractPlotModule(settings)
 437 {
 438     for (int i = 0; i < DIM; ++i)
 439     {
 440         bWrite_[i] = true;
 441     }
 442     bWrite_[DIM] = false;
 443 }
 444
 445
 446 void
 447 AnalysisDataVectorPlotModule::setWriteX(bool bWrite)
 448 {
 449     bWrite_[XX] = bWrite;
 450 }
 451
 452
 453 void
 454 AnalysisDataVectorPlotModule::setWriteY(bool bWrite)
 455 {
 456     bWrite_[YY] = bWrite;
 457 }
 458
 459
 460 void
 461 AnalysisDataVectorPlotModule::setWriteZ(bool bWrite)
 462 {
 463     bWrite_[ZZ] = bWrite;
 464 }
 465
 466
 467 void
 468 AnalysisDataVectorPlotModule::setWriteNorm(bool bWrite)
 469 {
 470     bWrite_[DIM] = bWrite;
 471 }
 472
 473
 474 void
 475 AnalysisDataVectorPlotModule::setWriteMask(bool bWrite[DIM + 1])
 476 {
 477     for (int i = 0; i < DIM + 1; ++i)
 478     {
 479         bWrite_[i] = bWrite[i];
 480     }
 481 }
 482
 483
 484 void
 485 AnalysisDataVectorPlotModule::pointsAdded(const AnalysisDataPointSetRef &points)
 486 {
 487     if (points.firstColumn() % DIM != 0)
 488     {
 489         GMX_THROW(APIError("Partial data points"));
 490     }
 491     if (!isFileOpen())
 492     {
 493         return;
 494     }
 495     for (int i = 0; i < points.columnCount(); i += 3)
 496     {
 497         for (int d = 0; d < DIM; ++d)
 498         {
 499             if (bWrite_[i])
 500             {
 501                 writeValue(points.y(i + d));
 502             }
 503         }
 504         if (bWrite_[DIM])
 505         {
 506             rvec y = { points.y(i), points.y(i + 1), points.y(i + 2) };
 507             writeValue(norm(y));
 508         }
 509     }
 510 }
 511
 512 } // namespace gmx