Bulk change to const ref for mtop
[alexxy/gromacs.git] / src / gromacs / selection / sm_simple.cpp
1 /*
2  * This file is part of the GROMACS molecular simulation package.
3  *
4  * Copyright (c) 2009,2010,2011,2012,2013 by the GROMACS development team.
5  * Copyright (c) 2014,2015,2016,2017,2018 by the GROMACS development team.
6  * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
7  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8  * and including many others, as listed in the AUTHORS file in the
9  * top-level source directory and at http://www.gromacs.org.
10  *
11  * GROMACS is free software; you can redistribute it and/or
12  * modify it under the terms of the GNU Lesser General Public License
13  * as published by the Free Software Foundation; either version 2.1
14  * of the License, or (at your option) any later version.
15  *
16  * GROMACS is distributed in the hope that it will be useful,
17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19  * Lesser General Public License for more details.
20  *
21  * You should have received a copy of the GNU Lesser General Public
22  * License along with GROMACS; if not, see
23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
25  *
26  * If you want to redistribute modifications to GROMACS, please
27  * consider that scientific software is very special. Version
28  * control is crucial - bugs must be traceable. We will be happy to
29  * consider code for inclusion in the official distribution, but
30  * derived work must not be called official GROMACS. Details are found
31  * in the README & COPYING files - if they are missing, get the
32  * official version at http://www.gromacs.org.
33  *
34  * To help us fund GROMACS development, we humbly ask that you cite
35  * the research papers on the package. Check out http://www.gromacs.org.
36  */
37 /*! \internal \file
38  * \brief
39  * Implements simple keyword selection methods.
40  *
41  * \author Teemu Murtola <teemu.murtola@gmail.com>
42  * \ingroup module_selection
43  */
44 #include "gmxpre.h"
45
46 #include <cctype>
47
48 #include "gromacs/topology/mtop_lookup.h"
49 #include "gromacs/topology/topology.h"
50 #include "gromacs/utility/arraysize.h"
51 #include "gromacs/utility/exceptions.h"
52 #include "gromacs/utility/gmxassert.h"
53
54 #include "position.h"
55 #include "selmethod.h"
56 #include "selmethod_impl.h"
57
58 /** Evaluates the \p all selection keyword. */
59 static void evaluate_all(const gmx::SelMethodEvalContext& context,
60                          gmx_ana_index_t*                 g,
61                          gmx_ana_selvalue_t*              out,
62                          void*                            data);
63 /** Evaluates the \p none selection keyword. */
64 static void evaluate_none(const gmx::SelMethodEvalContext& context,
65                           gmx_ana_index_t*                 g,
66                           gmx_ana_selvalue_t*              out,
67                           void*                            data);
68 /** Evaluates the \p atomnr selection keyword. */
69 static void evaluate_atomnr(const gmx::SelMethodEvalContext& context,
70                             gmx_ana_index_t*                 g,
71                             gmx_ana_selvalue_t*              out,
72                             void*                            data);
73 /** Evaluates the \p resnr selection keyword. */
74 static void evaluate_resnr(const gmx::SelMethodEvalContext& context,
75                            gmx_ana_index_t*                 g,
76                            gmx_ana_selvalue_t*              out,
77                            void*                            data);
78 /** Evaluates the \p resindex selection keyword. */
79 static void evaluate_resindex(const gmx::SelMethodEvalContext& context,
80                               gmx_ana_index_t*                 g,
81                               gmx_ana_selvalue_t*              out,
82                               void*                            data);
83 /*! \brief
84  * Checks whether molecule information is present in the topology.
85  *
86  * \param[in] top  Topology structure.
87  * \param     npar Not used.
88  * \param     param Not used.
89  * \param     data Not used.
90  * \returns   0 if molecule info is present in the topology, -1 otherwise.
91  *
92  * If molecule information is not found, also prints an error message.
93  */
94 static void check_molecules(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
95 /** Evaluates the \p molindex selection keyword. */
96 static void evaluate_molindex(const gmx::SelMethodEvalContext& context,
97                               gmx_ana_index_t*                 g,
98                               gmx_ana_selvalue_t*              out,
99                               void*                            data);
100 /** Evaluates the \p atomname selection keyword. */
101 static void evaluate_atomname(const gmx::SelMethodEvalContext& context,
102                               gmx_ana_index_t*                 g,
103                               gmx_ana_selvalue_t*              out,
104                               void*                            data);
105 /** Evaluates the \p pdbatomname selection keyword. */
106 static void evaluate_pdbatomname(const gmx::SelMethodEvalContext& context,
107                                  gmx_ana_index_t*                 g,
108                                  gmx_ana_selvalue_t*              out,
109                                  void*                            data);
110 /*! \brief
111  * Checks whether atom types are present in the topology.
112  *
113  * \param[in] top  Topology structure.
114  * \param     npar Not used.
115  * \param     param Not used.
116  * \param     data Not used.
117  */
118 static void check_atomtype(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
119 /** Evaluates the \p atomtype selection keyword. */
120 static void evaluate_atomtype(const gmx::SelMethodEvalContext& context,
121                               gmx_ana_index_t*                 g,
122                               gmx_ana_selvalue_t*              out,
123                               void*                            data);
124 /** Evaluates the \p insertcode selection keyword. */
125 static void evaluate_insertcode(const gmx::SelMethodEvalContext& context,
126                                 gmx_ana_index_t*                 g,
127                                 gmx_ana_selvalue_t*              out,
128                                 void*                            data);
129 /** Evaluates the \p chain selection keyword. */
130 static void evaluate_chain(const gmx::SelMethodEvalContext& context,
131                            gmx_ana_index_t*                 g,
132                            gmx_ana_selvalue_t*              out,
133                            void*                            data);
134 /** Evaluates the \p mass selection keyword. */
135 static void evaluate_mass(const gmx::SelMethodEvalContext& context,
136                           gmx_ana_index_t*                 g,
137                           gmx_ana_selvalue_t*              out,
138                           void*                            data);
139 /*! \brief
140  * Checks whether charges are present in the topology.
141  *
142  * \param[in] top  Topology structure.
143  * \param     npar Not used.
144  * \param     param Not used.
145  * \param     data Not used.
146  */
147 static void check_charge(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
148 /** Evaluates the \p charge selection keyword. */
149 static void evaluate_charge(const gmx::SelMethodEvalContext& context,
150                             gmx_ana_index_t*                 g,
151                             gmx_ana_selvalue_t*              out,
152                             void*                            data);
153 /*! \brief
154  * Checks whether PDB info is present in the topology.
155  *
156  * \param[in] top  Topology structure.
157  * \param     npar Not used.
158  * \param     param Not used.
159  * \param     data Not used.
160  * \returns   0 if PDB info is present in the topology, -1 otherwise.
161  *
162  * If PDB info is not found, also prints an error message.
163  */
164 static void check_pdbinfo(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
165 /** Evaluates the \p altloc selection keyword. */
166 static void evaluate_altloc(const gmx::SelMethodEvalContext& context,
167                             gmx_ana_index_t*                 g,
168                             gmx_ana_selvalue_t*              out,
169                             void*                            data);
170 /** Evaluates the \p occupancy selection keyword. */
171 static void evaluate_occupancy(const gmx::SelMethodEvalContext& context,
172                                gmx_ana_index_t*                 g,
173                                gmx_ana_selvalue_t*              out,
174                                void*                            data);
175 /** Evaluates the \p betafactor selection keyword. */
176 static void evaluate_betafactor(const gmx::SelMethodEvalContext& context,
177                                 gmx_ana_index_t*                 g,
178                                 gmx_ana_selvalue_t*              out,
179                                 void*                            data);
180 /** Evaluates the \p resname selection keyword. */
181 static void evaluate_resname(const gmx::SelMethodEvalContext& context,
182                              gmx_ana_index_t*                 g,
183                              gmx_ana_selvalue_t*              out,
184                              void*                            data);
185
186 /** Evaluates the \p x selection keyword. */
187 static void evaluate_x(const gmx::SelMethodEvalContext& context,
188                        gmx_ana_pos_t*                   pos,
189                        gmx_ana_selvalue_t*              out,
190                        void*                            data);
191 /** Evaluates the \p y selection keyword. */
192 static void evaluate_y(const gmx::SelMethodEvalContext& context,
193                        gmx_ana_pos_t*                   pos,
194                        gmx_ana_selvalue_t*              out,
195                        void*                            data);
196 /** Evaluates the \p z selection keyword. */
197 static void evaluate_z(const gmx::SelMethodEvalContext& context,
198                        gmx_ana_pos_t*                   pos,
199                        gmx_ana_selvalue_t*              out,
200                        void*                            data);
201
202 //! Help title for atom name selection keywords.
203 static const char helptitle_atomname[] = "Selecting atoms by name";
204 //! Help text for atom name selection keywords.
205 static const char* const help_atomname[] = {
206     "::",
207     "",
208     "  name",
209     "  pdbname",
210     "  atomname",
211     "  pdbatomname",
212     "",
213     "These keywords select atoms by name. [TT]name[tt] selects atoms using",
214     "the GROMACS atom naming convention.",
215     "For input formats other than PDB, the atom names are matched exactly",
216     "as they appear in the input file. For PDB files, 4 character atom names",
217     "that start with a digit are matched after moving the digit to the end",
218     "(e.g., to match 3HG2 from a PDB file, use [TT]name HG23[tt]).",
219     "[TT]pdbname[tt] can only be used with a PDB input file, and selects",
220     "atoms based on the exact name given in the input file, without the",
221     "transformation described above.[PAR]",
222
223     "[TT]atomname[tt] and [TT]pdbatomname[tt] are synonyms for the above two",
224     "keywords."
225 };
226
227 //! Help title for residue index selection keywords.
228 static const char helptitle_resindex[] = "Selecting atoms by residue number";
229 //! Help text for residue index selection keywords.
230 static const char* const help_resindex[] = {
231     "::",
232     "",
233     "  resnr",
234     "  resid",
235     "  resindex",
236     "  residue",
237     "",
238     "[TT]resnr[tt] selects atoms using the residue numbering in the input",
239     "file. [TT]resid[tt] is synonym for this keyword for VMD compatibility.",
240     "",
241     "[TT]resindex N[tt] selects the [TT]N[tt] th residue starting from the",
242     "beginning of the input file. This is useful for uniquely identifying",
243     "residues if there are duplicate numbers in the input file (e.g., in",
244     "multiple chains).",
245     "[TT]residue[tt] is a synonym for [TT]resindex[tt]. This allows",
246     "[TT]same residue as[tt] to work as expected."
247 };
248
249 /** Selection method data for \p all selection keyword. */
250 gmx_ana_selmethod_t sm_all = {
251     "all",   GROUP_VALUE, 0,       0,       nullptr,       nullptr, nullptr,
252     nullptr, nullptr,     nullptr, nullptr, &evaluate_all, nullptr,
253 };
254
255 /** Selection method data for \p none selection keyword. */
256 gmx_ana_selmethod_t sm_none = {
257     "none",  GROUP_VALUE, 0,       0,       nullptr,        nullptr, nullptr,
258     nullptr, nullptr,     nullptr, nullptr, &evaluate_none, nullptr,
259 };
260
261 /** Selection method data for \p atomnr selection keyword. */
262 gmx_ana_selmethod_t sm_atomnr = {
263     "atomnr", INT_VALUE, 0,       0,       nullptr,          nullptr, nullptr,
264     nullptr,  nullptr,   nullptr, nullptr, &evaluate_atomnr, nullptr,
265 };
266
267 /** Selection method data for \p resnr selection keyword. */
268 gmx_ana_selmethod_t sm_resnr = {
269     "resnr",      INT_VALUE,
270     SMETH_REQTOP, 0,
271     nullptr,      nullptr,
272     nullptr,      nullptr,
273     nullptr,      nullptr,
274     nullptr,      &evaluate_resnr,
275     nullptr,      { nullptr, helptitle_resindex, asize(help_resindex), help_resindex }
276 };
277
278 /** Selection method data for \p resindex selection keyword. */
279 gmx_ana_selmethod_t sm_resindex = {
280     "resindex",   INT_VALUE,
281     SMETH_REQTOP, 0,
282     nullptr,      nullptr,
283     nullptr,      nullptr,
284     nullptr,      nullptr,
285     nullptr,      &evaluate_resindex,
286     nullptr,      { nullptr, helptitle_resindex, asize(help_resindex), help_resindex }
287 };
288
289 /** Selection method data for \p molindex selection keyword. */
290 gmx_ana_selmethod_t sm_molindex = {
291     "molindex", INT_VALUE,        SMETH_REQTOP, 0,       nullptr, nullptr,
292     nullptr,    &check_molecules, nullptr,      nullptr, nullptr, &evaluate_molindex,
293     nullptr,
294 };
295
296 /** Selection method data for \p atomname selection keyword. */
297 gmx_ana_selmethod_t sm_atomname = {
298     "atomname",   STR_VALUE,
299     SMETH_REQTOP, 0,
300     nullptr,      nullptr,
301     nullptr,      nullptr,
302     nullptr,      nullptr,
303     nullptr,      &evaluate_atomname,
304     nullptr,      { nullptr, helptitle_atomname, asize(help_atomname), help_atomname }
305 };
306
307 /** Selection method data for \p pdbatomname selection keyword. */
308 gmx_ana_selmethod_t sm_pdbatomname = {
309     "pdbatomname", STR_VALUE,
310     SMETH_REQTOP,  0,
311     nullptr,       nullptr,
312     nullptr,       &check_pdbinfo,
313     nullptr,       nullptr,
314     nullptr,       &evaluate_pdbatomname,
315     nullptr,       { nullptr, helptitle_atomname, asize(help_atomname), help_atomname }
316 };
317
318 /** Selection method data for \p atomtype selection keyword. */
319 gmx_ana_selmethod_t sm_atomtype = {
320     "atomtype", STR_VALUE,       SMETH_REQTOP, 0,       nullptr, nullptr,
321     nullptr,    &check_atomtype, nullptr,      nullptr, nullptr, &evaluate_atomtype,
322     nullptr,
323 };
324
325 /** Selection method data for \p resname selection keyword. */
326 gmx_ana_selmethod_t sm_resname = {
327     "resname", STR_VALUE, SMETH_REQTOP, 0,       nullptr,           nullptr, nullptr,
328     nullptr,   nullptr,   nullptr,      nullptr, &evaluate_resname, nullptr,
329 };
330
331 /** Selection method data for \p chain selection keyword. */
332 gmx_ana_selmethod_t sm_insertcode = {
333     "insertcode",
334     STR_VALUE,
335     SMETH_REQTOP | SMETH_CHARVAL,
336     0,
337     nullptr,
338     nullptr,
339     nullptr,
340     nullptr,
341     nullptr,
342     nullptr,
343     nullptr,
344     &evaluate_insertcode,
345     nullptr,
346 };
347
348 /** Selection method data for \p chain selection keyword. */
349 gmx_ana_selmethod_t sm_chain = {
350     "chain", STR_VALUE, SMETH_REQTOP | SMETH_CHARVAL,
351     0,       nullptr,   nullptr,
352     nullptr, nullptr,   nullptr,
353     nullptr, nullptr,   &evaluate_chain,
354     nullptr,
355 };
356
357 /** Selection method data for \p mass selection keyword. */
358 gmx_ana_selmethod_t sm_mass = {
359     "mass",  REAL_VALUE, SMETH_REQMASS, 0,       nullptr,        nullptr, nullptr,
360     nullptr, nullptr,    nullptr,       nullptr, &evaluate_mass, nullptr,
361 };
362
363 /** Selection method data for \p charge selection keyword. */
364 gmx_ana_selmethod_t sm_charge = {
365     "charge",      REAL_VALUE, SMETH_REQTOP, 0,       nullptr,          nullptr, nullptr,
366     &check_charge, nullptr,    nullptr,      nullptr, &evaluate_charge, nullptr,
367 };
368
369 /** Selection method data for \p chain selection keyword. */
370 gmx_ana_selmethod_t sm_altloc = {
371     "altloc", STR_VALUE,      SMETH_REQTOP | SMETH_CHARVAL,
372     0,        nullptr,        nullptr,
373     nullptr,  &check_pdbinfo, nullptr,
374     nullptr,  nullptr,        &evaluate_altloc,
375     nullptr,
376 };
377
378 /** Selection method data for \p occupancy selection keyword. */
379 gmx_ana_selmethod_t sm_occupancy = {
380     "occupancy", REAL_VALUE,     SMETH_REQTOP, 0,       nullptr, nullptr,
381     nullptr,     &check_pdbinfo, nullptr,      nullptr, nullptr, &evaluate_occupancy,
382     nullptr,
383 };
384
385 /** Selection method data for \p betafactor selection keyword. */
386 gmx_ana_selmethod_t sm_betafactor = {
387     "betafactor", REAL_VALUE,     SMETH_REQTOP, 0,       nullptr, nullptr,
388     nullptr,      &check_pdbinfo, nullptr,      nullptr, nullptr, &evaluate_betafactor,
389     nullptr,
390 };
391
392 /** Selection method data for \p x selection keyword. */
393 gmx_ana_selmethod_t sm_x = {
394     "x",     REAL_VALUE, SMETH_DYNAMIC, 0,       nullptr, nullptr,     nullptr,
395     nullptr, nullptr,    nullptr,       nullptr, nullptr, &evaluate_x,
396 };
397
398 /** Selection method data for \p y selection keyword. */
399 gmx_ana_selmethod_t sm_y = {
400     "y",     REAL_VALUE, SMETH_DYNAMIC, 0,       nullptr, nullptr,     nullptr,
401     nullptr, nullptr,    nullptr,       nullptr, nullptr, &evaluate_y,
402 };
403
404 /** Selection method data for \p z selection keyword. */
405 gmx_ana_selmethod_t sm_z = {
406     "z",     REAL_VALUE, SMETH_DYNAMIC, 0,       nullptr, nullptr,     nullptr,
407     nullptr, nullptr,    nullptr,       nullptr, nullptr, &evaluate_z,
408 };
409
410 /*!
411  * See sel_updatefunc() for description of the parameters.
412  * \p data is not used.
413  *
414  * Copies \p g to \p out->u.g.
415  */
416 static void evaluate_all(const gmx::SelMethodEvalContext& /*context*/,
417                          gmx_ana_index_t*    g,
418                          gmx_ana_selvalue_t* out,
419                          void* /* data */)
420 {
421     gmx_ana_index_copy(out->u.g, g, false);
422 }
423
424 /*!
425  * See sel_updatefunc() for description of the parameters.
426  * \p data is not used.
427  *
428  * Returns an empty \p out->u.g.
429  */
430 static void evaluate_none(const gmx::SelMethodEvalContext& /*context*/,
431                           gmx_ana_index_t* /* g */,
432                           gmx_ana_selvalue_t* out,
433                           void* /* data */)
434 {
435     out->u.g->isize = 0;
436 }
437
438 /*!
439  * See sel_updatefunc() for description of the parameters.
440  * \p data is not used.
441  *
442  * Returns the indices for each atom in \p out->u.i.
443  */
444 static void evaluate_atomnr(const gmx::SelMethodEvalContext& /*context*/,
445                             gmx_ana_index_t*    g,
446                             gmx_ana_selvalue_t* out,
447                             void* /* data */)
448 {
449     int i;
450
451     out->nr = g->isize;
452     for (i = 0; i < g->isize; ++i)
453     {
454         out->u.i[i] = g->index[i] + 1;
455     }
456 }
457
458 /*!
459  * See sel_updatefunc() for description of the parameters.
460  * \p data is not used.
461  *
462  * Returns the residue numbers for each atom in \p out->u.i.
463  */
464 static void evaluate_resnr(const gmx::SelMethodEvalContext& context,
465                            gmx_ana_index_t*                 g,
466                            gmx_ana_selvalue_t*              out,
467                            void* /* data */)
468 {
469     out->nr  = g->isize;
470     int molb = 0;
471     for (int i = 0; i < g->isize; ++i)
472     {
473         mtopGetAtomAndResidueName(*context.top, g->index[i], &molb, nullptr, &out->u.i[i], nullptr, nullptr);
474     }
475 }
476
477 /*!
478  * See sel_updatefunc() for description of the parameters.
479  * \p data is not used.
480  *
481  * Returns the residue indices for each atom in \p out->u.i.
482  */
483 static void evaluate_resindex(const gmx::SelMethodEvalContext& context,
484                               gmx_ana_index_t*                 g,
485                               gmx_ana_selvalue_t*              out,
486                               void* /* data */)
487 {
488     out->nr  = g->isize;
489     int molb = 0;
490     for (int i = 0; i < g->isize; ++i)
491     {
492         int resind;
493         mtopGetAtomAndResidueName(*context.top, g->index[i], &molb, nullptr, nullptr, nullptr, &resind);
494         out->u.i[i] = resind + 1;
495     }
496 }
497
498 static void check_molecules(const gmx_mtop_t* top, int /* npar */, gmx_ana_selparam_t* /* param */, void* /* data */)
499 {
500     bool bOk;
501
502     bOk = (top != nullptr && top->haveMoleculeIndices);
503     if (!bOk)
504     {
505         GMX_THROW(gmx::InconsistentInputError("Molecule information not available in topology"));
506     }
507 }
508
509 /*!
510  * See sel_updatefunc() for description of the parameters.
511  * \p data is not used.
512  *
513  * Returns the molecule indices for each atom in \p out->u.i.
514  */
515 static void evaluate_molindex(const gmx::SelMethodEvalContext& context,
516                               gmx_ana_index_t*                 g,
517                               gmx_ana_selvalue_t*              out,
518                               void* /* data */)
519 {
520     out->nr  = g->isize;
521     int molb = 0;
522     for (int i = 0; i < g->isize; ++i)
523     {
524         out->u.i[i] = mtopGetMoleculeIndex(*context.top, g->index[i], &molb) + 1;
525     }
526 }
527
528 /*!
529  * See sel_updatefunc() for description of the parameters.
530  * \p data is not used.
531  *
532  * Returns the atom name for each atom in \p out->u.s.
533  */
534 static void evaluate_atomname(const gmx::SelMethodEvalContext& context,
535                               gmx_ana_index_t*                 g,
536                               gmx_ana_selvalue_t*              out,
537                               void* /* data */)
538 {
539     out->nr  = g->isize;
540     int molb = 0;
541     for (int i = 0; i < g->isize; ++i)
542     {
543         const char* atom_name;
544         mtopGetAtomAndResidueName(*context.top, g->index[i], &molb, &atom_name, nullptr, nullptr, nullptr);
545         out->u.s[i] = const_cast<char*>(atom_name);
546     }
547 }
548
549 /*!
550  * See sel_updatefunc() for description of the parameters.
551  * \p data is not used.
552  *
553  * Returns the PDB atom name for each atom in \p out->u.s.
554  */
555 static void evaluate_pdbatomname(const gmx::SelMethodEvalContext& context,
556                                  gmx_ana_index_t*                 g,
557                                  gmx_ana_selvalue_t*              out,
558                                  void* /* data */)
559 {
560     out->nr  = g->isize;
561     int molb = 0;
562     for (int i = 0; i < g->isize; ++i)
563     {
564         const char* s = mtopGetAtomPdbInfo(*context.top, g->index[i], &molb).atomnm;
565         while (std::isspace(*s))
566         {
567             ++s;
568         }
569         out->u.s[i] = const_cast<char*>(s);
570     }
571 }
572
573 static void check_atomtype(const gmx_mtop_t* top, int /* npar */, gmx_ana_selparam_t* /* param */, void* /* data */)
574 {
575     if (!gmx_mtop_has_atomtypes(top))
576     {
577         GMX_THROW(gmx::InconsistentInputError("Atom types not available in topology"));
578     }
579 }
580
581 /*!
582  * See sel_updatefunc() for description of the parameters.
583  * \p data is not used.
584  *
585  * Returns the atom type for each atom in \p out->u.s.
586  * Segfaults if atom types are not found in the topology.
587  */
588 static void evaluate_atomtype(const gmx::SelMethodEvalContext& context,
589                               gmx_ana_index_t*                 g,
590                               gmx_ana_selvalue_t*              out,
591                               void* /* data */)
592 {
593     out->nr  = g->isize;
594     int molb = 0;
595     for (int i = 0; i < g->isize; ++i)
596     {
597         int atomIndexInMolecule;
598         mtopGetMolblockIndex(*context.top, g->index[i], &molb, nullptr, &atomIndexInMolecule);
599         const gmx_moltype_t& moltype = context.top->moltype[context.top->molblock[molb].type];
600         out->u.s[i]                  = *moltype.atoms.atomtype[atomIndexInMolecule];
601     }
602 }
603
604 /*!
605  * See sel_updatefunc() for description of the parameters.
606  * \p data is not used.
607  *
608  * Returns the residue name for each atom in \p out->u.s.
609  */
610 static void evaluate_resname(const gmx::SelMethodEvalContext& context,
611                              gmx_ana_index_t*                 g,
612                              gmx_ana_selvalue_t*              out,
613                              void* /* data */)
614 {
615     out->nr  = g->isize;
616     int molb = 0;
617     for (int i = 0; i < g->isize; ++i)
618     {
619         out->u.s[i] = *mtopGetResidueInfo(*context.top, g->index[i], &molb).name;
620     }
621 }
622
623 /*!
624  * See sel_updatefunc() for description of the parameters.
625  * \p data is not used.
626  *
627  * Returns the insertion code for each atom in \p out->u.s.
628  */
629 static void evaluate_insertcode(const gmx::SelMethodEvalContext& context,
630                                 gmx_ana_index_t*                 g,
631                                 gmx_ana_selvalue_t*              out,
632                                 void* /* data */)
633 {
634     out->nr  = g->isize;
635     int molb = 0;
636     for (int i = 0; i < g->isize; ++i)
637     {
638         out->u.s[i][0] = mtopGetResidueInfo(*context.top, g->index[i], &molb).ic;
639     }
640 }
641
642 /*!
643  * See sel_updatefunc() for description of the parameters.
644  * \p data is not used.
645  *
646  * Returns the chain for each atom in \p out->u.s.
647  */
648 static void evaluate_chain(const gmx::SelMethodEvalContext& context,
649                            gmx_ana_index_t*                 g,
650                            gmx_ana_selvalue_t*              out,
651                            void* /* data */)
652 {
653     out->nr  = g->isize;
654     int molb = 0;
655     for (int i = 0; i < g->isize; ++i)
656     {
657         out->u.s[i][0] = mtopGetResidueInfo(*context.top, g->index[i], &molb).chainid;
658     }
659 }
660
661 /*!
662  * See sel_updatefunc() for description of the parameters.
663  * \p data is not used.
664  *
665  * Returns the mass for each atom in \p out->u.r.
666  */
667 static void evaluate_mass(const gmx::SelMethodEvalContext& context,
668                           gmx_ana_index_t*                 g,
669                           gmx_ana_selvalue_t*              out,
670                           void* /* data */)
671 {
672     GMX_RELEASE_ASSERT(gmx_mtop_has_masses(context.top), "Masses not available for evaluation");
673     out->nr  = g->isize;
674     int molb = 0;
675     for (int i = 0; i < g->isize; ++i)
676     {
677         out->u.r[i] = mtopGetAtomMass(*context.top, g->index[i], &molb);
678     }
679 }
680
681
682 static void check_charge(const gmx_mtop_t* top, int /* npar */, gmx_ana_selparam_t* /* param */, void* /* data */)
683 {
684     if (!gmx_mtop_has_charges(top))
685     {
686         GMX_THROW(gmx::InconsistentInputError("Charges not available in topology"));
687     }
688 }
689
690 /*!
691  * See sel_updatefunc() for description of the parameters.
692  * \p data is not used.
693  *
694  * Returns the charge for each atom in \p out->u.r.
695  */
696 static void evaluate_charge(const gmx::SelMethodEvalContext& context,
697                             gmx_ana_index_t*                 g,
698                             gmx_ana_selvalue_t*              out,
699                             void* /* data */)
700 {
701     out->nr  = g->isize;
702     int molb = 0;
703     for (int i = 0; i < g->isize; ++i)
704     {
705         out->u.r[i] = mtopGetAtomParameters(*context.top, g->index[i], &molb).q;
706     }
707 }
708
709 static void check_pdbinfo(const gmx_mtop_t* top, int /* npar */, gmx_ana_selparam_t* /* param */, void* /* data */)
710 {
711     if (!gmx_mtop_has_pdbinfo(top))
712     {
713         GMX_THROW(gmx::InconsistentInputError("PDB info not available in topology"));
714     }
715 }
716
717 /*!
718  * See sel_updatefunc() for description of the parameters.
719  * \p data is not used.
720  *
721  * Returns the alternate location identifier for each atom in \p out->u.s.
722  */
723 static void evaluate_altloc(const gmx::SelMethodEvalContext& context,
724                             gmx_ana_index_t*                 g,
725                             gmx_ana_selvalue_t*              out,
726                             void* /* data */)
727 {
728     out->nr  = g->isize;
729     int molb = 0;
730     for (int i = 0; i < g->isize; ++i)
731     {
732         out->u.s[i][0] = mtopGetAtomPdbInfo(*context.top, g->index[i], &molb).altloc;
733     }
734 }
735
736 /*!
737  * See sel_updatefunc() for description of the parameters.
738  * \p data is not used.
739  *
740  * Returns the occupancy numbers for each atom in \p out->u.r.
741  * Segfaults if PDB info is not found in the topology.
742  */
743 static void evaluate_occupancy(const gmx::SelMethodEvalContext& context,
744                                gmx_ana_index_t*                 g,
745                                gmx_ana_selvalue_t*              out,
746                                void* /* data */)
747 {
748     out->nr  = g->isize;
749     int molb = 0;
750     for (int i = 0; i < g->isize; ++i)
751     {
752         out->u.r[i] = mtopGetAtomPdbInfo(*context.top, g->index[i], &molb).occup;
753     }
754 }
755
756 /*!
757  * See sel_updatefunc() for description of the parameters.
758  * \p data is not used.
759  *
760  * Returns the B-factors for each atom in \p out->u.r.
761  * Segfaults if PDB info is not found in the topology.
762  */
763 static void evaluate_betafactor(const gmx::SelMethodEvalContext& context,
764                                 gmx_ana_index_t*                 g,
765                                 gmx_ana_selvalue_t*              out,
766                                 void* /* data */)
767 {
768     out->nr  = g->isize;
769     int molb = 0;
770     for (int i = 0; i < g->isize; ++i)
771     {
772         out->u.r[i] = mtopGetAtomPdbInfo(*context.top, g->index[i], &molb).bfac;
773     }
774 }
775
776 /*! \brief
777  * Internal utility function for position keyword evaluation.
778  *
779  * \param[out] out  Output array.
780  * \param[in]  pos  Position data to use instead of atomic coordinates.
781  * \param[in]  d    Coordinate index to evaluate (\p XX, \p YY or \p ZZ).
782  *
783  * This function is used internally by evaluate_x(), evaluate_y() and
784  * evaluate_z() to do the actual evaluation.
785  */
786 static void evaluate_coord(real out[], gmx_ana_pos_t* pos, int d)
787 {
788     for (int i = 0; i < pos->count(); ++i)
789     {
790         out[i] = pos->x[i][d];
791     }
792     // TODO: Make this more efficient by directly extracting the coordinates
793     // from the frame coordinates for atomic positions instead of going through
794     // a position calculation.
795 }
796
797 /*!
798  * See sel_updatefunc_pos() for description of the parameters.
799  * \p data is not used.
800  *
801  * Returns the \p x coordinate for each position in \p out->u.r.
802  */
803 static void evaluate_x(const gmx::SelMethodEvalContext& /*context*/,
804                        gmx_ana_pos_t*      pos,
805                        gmx_ana_selvalue_t* out,
806                        void* /*data*/)
807 {
808     out->nr = pos->count();
809     evaluate_coord(out->u.r, pos, XX);
810 }
811
812 /*!
813  * See sel_updatefunc() for description of the parameters.
814  * \p data is not used.
815  *
816  * Returns the \p y coordinate for each position in \p out->u.r.
817  */
818 static void evaluate_y(const gmx::SelMethodEvalContext& /*context*/,
819                        gmx_ana_pos_t*      pos,
820                        gmx_ana_selvalue_t* out,
821                        void* /*data*/)
822 {
823     out->nr = pos->count();
824     evaluate_coord(out->u.r, pos, YY);
825 }
826
827 /*!
828  * See sel_updatefunc() for description of the parameters.
829  * \p data is not used.
830  *
831  * Returns the \p z coordinate for each position in \p out->u.r.
832  */
833 static void evaluate_z(const gmx::SelMethodEvalContext& /*context*/,
834                        gmx_ana_pos_t*      pos,
835                        gmx_ana_selvalue_t* out,
836                        void* /*data*/)
837 {
838     out->nr = pos->count();
839     evaluate_coord(out->u.r, pos, ZZ);
840 }