2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2010-2018, The GROMACS development team.
5 * Copyright (c) 2019,2020, by the GROMACS development team, led by
6 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
7 * and including many others, as listed in the AUTHORS file in the
8 * top-level source directory and at http://www.gromacs.org.
10 * GROMACS is free software; you can redistribute it and/or
11 * modify it under the terms of the GNU Lesser General Public License
12 * as published by the Free Software Foundation; either version 2.1
13 * of the License, or (at your option) any later version.
15 * GROMACS is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 * Lesser General Public License for more details.
20 * You should have received a copy of the GNU Lesser General Public
21 * License along with GROMACS; if not, see
22 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
23 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
25 * If you want to redistribute modifications to GROMACS, please
26 * consider that scientific software is very special. Version
27 * control is crucial - bugs must be traceable. We will be happy to
28 * consider code for inclusion in the official distribution, but
29 * derived work must not be called official GROMACS. Details are found
30 * in the README & COPYING files - if they are missing, get the
31 * official version at http://www.gromacs.org.
33 * To help us fund GROMACS development, we humbly ask that you cite
34 * the research papers on the package. Check out http://www.gromacs.org.
42 #include "gromacs/commandline/pargs.h"
43 #include "gromacs/fileio/confio.h"
44 #include "gromacs/fileio/matio.h"
45 #include "gromacs/fileio/trxio.h"
46 #include "gromacs/fileio/xvgr.h"
47 #include "gromacs/gmxana/binsearch.h"
48 #include "gromacs/gmxana/gmx_ana.h"
49 #include "gromacs/gmxana/gstat.h"
50 #include "gromacs/gmxana/powerspect.h"
51 #include "gromacs/math/vec.h"
52 #include "gromacs/pbcutil/pbc.h"
53 #include "gromacs/pbcutil/rmpbc.h"
54 #include "gromacs/topology/index.h"
55 #include "gromacs/topology/topology.h"
56 #include "gromacs/utility/arraysize.h"
57 #include "gromacs/utility/cstringutil.h"
58 #include "gromacs/utility/fatalerror.h"
59 #include "gromacs/utility/futil.h"
60 #include "gromacs/utility/gmxassert.h"
61 #include "gromacs/utility/smalloc.h"
63 static void find_tetra_order_grid(t_topology top
,
78 int ix
, jx
, i
, j
, k
, l
, n
, *nn
[4];
79 rvec dx
, rj
, rk
, urk
, urj
;
80 real cost
, cost2
, *sgmol
, *skmol
, rmean
, rmean2
, r2
, box2
, *r_nn
[4];
82 int slindex_x
, slindex_y
, slindex_z
;
84 real onethird
= 1.0 / 3.0;
87 /* dmat = init_mat(maxidx, FALSE); */
89 box2
= box
[XX
][XX
] * box
[XX
][XX
];
91 /* Initialize expanded sl_count array */
92 snew(sl_count
, nslicex
);
93 for (i
= 0; i
< nslicex
; i
++)
95 snew(sl_count
[i
], nslicey
);
96 for (j
= 0; j
< nslicey
; j
++)
98 snew(sl_count
[i
][j
], nslicez
);
103 for (i
= 0; (i
< 4); i
++)
105 snew(r_nn
[i
], natoms
);
108 for (j
= 0; (j
< natoms
); j
++)
117 /* Must init pbc every step because of pressure coupling */
118 set_pbc(&pbc
, pbcType
, box
);
119 gpbc
= gmx_rmpbc_init(&top
.idef
, pbcType
, natoms
);
120 gmx_rmpbc(gpbc
, natoms
, box
, x
);
125 for (i
= 0; (i
< maxidx
); i
++)
126 { /* loop over index file */
128 for (j
= 0; (j
< maxidx
); j
++)
138 pbc_dx(&pbc
, x
[ix
], x
[jx
], dx
);
141 /* set_mat_entry(dmat,i,j,r2); */
143 /* determine the nearest neighbours */
146 r_nn
[3][i
] = r_nn
[2][i
];
148 r_nn
[2][i
] = r_nn
[1][i
];
150 r_nn
[1][i
] = r_nn
[0][i
];
155 else if (r2
< r_nn
[1][i
])
157 r_nn
[3][i
] = r_nn
[2][i
];
159 r_nn
[2][i
] = r_nn
[1][i
];
164 else if (r2
< r_nn
[2][i
])
166 r_nn
[3][i
] = r_nn
[2][i
];
171 else if (r2
< r_nn
[3][i
])
179 /* calculate mean distance between nearest neighbours */
181 for (j
= 0; (j
< 4); j
++)
183 r_nn
[j
][i
] = std::sqrt(r_nn
[j
][i
]);
192 /* Chau1998a eqn 3 */
193 /* angular part tetrahedrality order parameter per atom */
194 for (j
= 0; (j
< 3); j
++)
196 for (k
= j
+ 1; (k
< 4); k
++)
198 pbc_dx(&pbc
, x
[ix
], x
[index
[nn
[k
][i
]]], rk
);
199 pbc_dx(&pbc
, x
[ix
], x
[index
[nn
[j
][i
]]], rj
);
204 cost
= iprod(urk
, urj
) + onethird
;
212 /* normalize sgmol between 0.0 and 1.0 */
213 sgmol
[i
] = 3 * sgmol
[i
] / 32;
216 /* distance part tetrahedrality order parameter per atom */
217 rmean2
= 4 * 3 * rmean
* rmean
;
218 for (j
= 0; (j
< 4); j
++)
220 skmol
[i
] += (rmean
- r_nn
[j
][i
]) * (rmean
- r_nn
[j
][i
]) / rmean2
;
221 /* printf("%d %f (%f %f %f %f) \n",
222 i, skmol[i], rmean, rmean2, r_nn[j][i], (rmean - r_nn[j][i]) );
228 /* Compute sliced stuff in x y z*/
229 slindex_x
= static_cast<int>(std::round((1 + x
[i
][XX
] / box
[XX
][XX
]) * nslicex
)) % nslicex
;
230 slindex_y
= static_cast<int>(std::round((1 + x
[i
][YY
] / box
[YY
][YY
]) * nslicey
)) % nslicey
;
231 slindex_z
= static_cast<int>(std::round((1 + x
[i
][ZZ
] / box
[ZZ
][ZZ
]) * nslicez
)) % nslicez
;
232 sggrid
[slindex_x
][slindex_y
][slindex_z
] += sgmol
[i
];
233 skgrid
[slindex_x
][slindex_y
][slindex_z
] += skmol
[i
];
234 (sl_count
[slindex_x
][slindex_y
][slindex_z
])++;
235 } /* loop over entries in index file */
240 for (i
= 0; (i
< nslicex
); i
++)
242 for (j
= 0; j
< nslicey
; j
++)
244 for (k
= 0; k
< nslicez
; k
++)
246 if (sl_count
[i
][j
][k
] > 0)
248 sggrid
[i
][j
][k
] /= sl_count
[i
][j
][k
];
249 skgrid
[i
][j
][k
] /= sl_count
[i
][j
][k
];
258 for (i
= 0; (i
< 4); i
++)
265 /*Determines interface from tetrahedral order parameter in box with specified binwidth. */
266 /*Outputs interface positions(bins), the number of timeframes, and the number of surface-mesh points in xy*/
268 static void calc_tetra_order_interface(const char* fnNDX
,
279 gmx_output_env_t
* oenv
)
281 FILE * fpsg
= nullptr, *fpsk
= nullptr;
290 int** index
= nullptr;
291 char** grpname
= nullptr;
292 int i
, j
, k
, n
, *isize
, ng
, nslicez
, framenr
;
293 real
***sg_grid
= nullptr, ***sk_grid
= nullptr, ***sg_fravg
= nullptr, ***sk_fravg
= nullptr,
294 ****sk_4d
= nullptr, ****sg_4d
= nullptr;
298 const real onehalf
= 1.0 / 2.0;
299 /* real ***intfpos[2]; pointers to arrays of two interface positions zcoord(framenr,xbin,ybin): intfpos[interface_index][t][nslicey*x+y]
300 * i.e 1D Row-major order in (t,x,y) */
303 read_tps_conf(fnTPS
, &top
, &pbcType
, &xtop
, nullptr, box
, FALSE
);
305 *nslicex
= static_cast<int>(box
[XX
][XX
] / binw
+ onehalf
); /*Calculate slicenr from binwidth*/
306 *nslicey
= static_cast<int>(box
[YY
][YY
] / binw
+ onehalf
);
307 nslicez
= static_cast<int>(box
[ZZ
][ZZ
] / binw
+ onehalf
);
311 /* get index groups */
312 printf("Select the group that contains the atoms you want to use for the tetrahedrality order "
313 "parameter calculation:\n");
317 get_index(&top
.atoms
, fnNDX
, ng
, isize
, index
, grpname
);
319 /* Analyze trajectory */
320 natoms
= read_first_x(oenv
, &status
, fnTRX
, &t
, &x
, box
);
321 if (natoms
> top
.atoms
.nr
)
323 gmx_fatal(FARGS
, "Topology (%d atoms) does not match trajectory (%d atoms)", top
.atoms
.nr
, natoms
);
325 check_index(nullptr, ng
, index
[0], nullptr, natoms
);
328 /*Prepare structures for temporary storage of frame info*/
329 snew(sg_grid
, *nslicex
);
330 snew(sk_grid
, *nslicex
);
331 for (i
= 0; i
< *nslicex
; i
++)
333 snew(sg_grid
[i
], *nslicey
);
334 snew(sk_grid
[i
], *nslicey
);
335 for (j
= 0; j
< *nslicey
; j
++)
337 snew(sg_grid
[i
][j
], nslicez
);
338 snew(sk_grid
[i
][j
], nslicez
);
347 /* Loop over frames*/
350 /*Initialize box meshes (temporary storage for each tblock frame -reinitialise every tblock steps */
351 if (framenr
% tblock
== 0)
353 srenew(sk_4d
, *nframes
+ 1);
354 srenew(sg_4d
, *nframes
+ 1);
355 snew(sg_fravg
, *nslicex
);
356 snew(sk_fravg
, *nslicex
);
357 for (i
= 0; i
< *nslicex
; i
++)
359 snew(sg_fravg
[i
], *nslicey
);
360 snew(sk_fravg
[i
], *nslicey
);
361 for (j
= 0; j
< *nslicey
; j
++)
363 snew(sg_fravg
[i
][j
], nslicez
);
364 snew(sk_fravg
[i
][j
], nslicez
);
369 find_tetra_order_grid(top
, pbcType
, natoms
, box
, x
, isize
[0], index
[0], &sg
, &sk
, *nslicex
,
370 *nslicey
, nslicez
, sg_grid
, sk_grid
);
371 GMX_RELEASE_ASSERT(sk_fravg
!= nullptr, "Trying to dereference NULL sk_fravg pointer");
372 for (i
= 0; i
< *nslicex
; i
++)
374 for (j
= 0; j
< *nslicey
; j
++)
376 for (k
= 0; k
< nslicez
; k
++)
378 sk_fravg
[i
][j
][k
] += sk_grid
[i
][j
][k
] / tblock
;
379 sg_fravg
[i
][j
][k
] += sg_grid
[i
][j
][k
] / tblock
;
386 if (framenr
% tblock
== 0)
388 GMX_RELEASE_ASSERT(sk_4d
!= nullptr, "Trying to dereference NULL sk_4d pointer");
389 sk_4d
[*nframes
] = sk_fravg
;
390 sg_4d
[*nframes
] = sg_fravg
;
394 } while (read_next_x(oenv
, status
, &t
, x
, box
));
401 /*Debugging for printing out the entire order parameter meshes.*/
404 fpsg
= xvgropen("sg_ang_mesh", "S\\sg\\N Angle Order Parameter / Meshpoint", "(nm)",
406 fpsk
= xvgropen("sk_dist_mesh", "S\\sk\\N Distance Order Parameter / Meshpoint", "(nm)",
408 for (n
= 0; n
< (*nframes
); n
++)
410 fprintf(fpsg
, "%i\n", n
);
411 fprintf(fpsk
, "%i\n", n
);
412 for (i
= 0; (i
< *nslicex
); i
++)
414 for (j
= 0; j
< *nslicey
; j
++)
416 for (k
= 0; k
< nslicez
; k
++)
418 fprintf(fpsg
, "%4f %4f %4f %8f\n", (i
+ 0.5) * box
[XX
][XX
] / (*nslicex
),
419 (j
+ 0.5) * box
[YY
][YY
] / (*nslicey
),
420 (k
+ 0.5) * box
[ZZ
][ZZ
] / nslicez
, sg_4d
[n
][i
][j
][k
]);
421 fprintf(fpsk
, "%4f %4f %4f %8f\n", (i
+ 0.5) * box
[XX
][XX
] / (*nslicex
),
422 (j
+ 0.5) * box
[YY
][YY
] / (*nslicey
),
423 (k
+ 0.5) * box
[ZZ
][ZZ
] / nslicez
, sk_4d
[n
][i
][j
][k
]);
433 /* Find positions of interface z by scanning orderparam for each frame and for each xy-mesh cylinder along z*/
435 /*Simple trial: assume interface is in the middle of -sgang1 and sgang2*/
436 sgintf
= 0.5 * (sgang1
+ sgang2
);
439 /*Allocate memory for interface arrays; */
441 snew((*intfpos
)[0], *nframes
);
442 snew((*intfpos
)[1], *nframes
);
444 bins
= (*nslicex
) * (*nslicey
);
447 snew(perm
, nslicez
); /*permutation array for sorting along normal coordinate*/
450 for (n
= 0; n
< *nframes
; n
++)
452 snew((*intfpos
)[0][n
], bins
);
453 snew((*intfpos
)[1][n
], bins
);
454 for (i
= 0; i
< *nslicex
; i
++)
456 for (j
= 0; j
< *nslicey
; j
++)
458 rangeArray(perm
, nslicez
); /*reset permutation array to identity*/
459 /*Binsearch returns 2 bin-numbers where the order param is <= setpoint sgintf*/
460 ndx1
= start_binsearch(sg_4d
[n
][i
][j
], perm
, 0, nslicez
/ 2 - 1, sgintf
, 1);
461 ndx2
= start_binsearch(sg_4d
[n
][i
][j
], perm
, nslicez
/ 2, nslicez
- 1, sgintf
, -1);
462 /*Use linear interpolation to smooth out the interface position*/
464 /*left interface (0)*/
465 /*if((sg_4d[n][i][j][perm[ndx1+1]]-sg_4d[n][i][j][perm[ndx1]])/sg_4d[n][i][j][perm[ndx1]] > 0.01){
466 pos=( (sgintf-sg_4d[n][i][j][perm[ndx1]])*perm[ndx1+1]+(sg_4d[n][i][j][perm[ndx1+1]]-sgintf)*perm[ndx1 ])*/
467 (*intfpos
)[0][n
][j
+ *nslicey
* i
] = (perm
[ndx1
] + onehalf
) * binw
;
468 /*right interface (1)*/
469 /*alpha=(sgintf-sg_4d[n][i][j][perm[ndx2]])/(sg_4d[n][i][j][perm[ndx2]+1]-sg_4d[n][i][j][perm[ndx2]]);*/
470 /*(*intfpos)[1][n][j+*nslicey*i]=((1-alpha)*perm[ndx2]+alpha*(perm[ndx2]+1)+onehalf)*box[ZZ][ZZ]/nslicez;*/
471 (*intfpos
)[1][n
][j
+ *nslicey
* i
] = (perm
[ndx2
] + onehalf
) * binw
;
481 static void writesurftoxpms(real
*** surf
,
486 gmx::ArrayRef
<const std::string
> outfiles
,
492 real
**profile1
, **profile2
;
493 real max1
, max2
, min1
, min2
, *xticks
, *yticks
;
494 t_rgb lo
= { 1, 1, 1 };
495 t_rgb hi
= { 0, 0, 0 };
496 FILE * xpmfile1
, *xpmfile2
;
498 /*Prepare xpm structures for output*/
500 /*Allocate memory to tick's and matrices*/
501 snew(xticks
, xbins
+ 1);
502 snew(yticks
, ybins
+ 1);
504 profile1
= mk_matrix(xbins
, ybins
, FALSE
);
505 profile2
= mk_matrix(xbins
, ybins
, FALSE
);
507 for (i
= 0; i
< xbins
+ 1; i
++)
511 for (j
= 0; j
< ybins
+ 1; j
++)
516 xpmfile1
= gmx_ffopen(outfiles
[0], "w");
517 xpmfile2
= gmx_ffopen(outfiles
[1], "w");
520 min1
= min2
= 1000.00;
522 for (n
= 0; n
< tblocks
; n
++)
524 sprintf(numbuf
, "%5d", n
);
525 /*Filling matrices for inclusion in xpm-files*/
526 for (i
= 0; i
< xbins
; i
++)
528 for (j
= 0; j
< ybins
; j
++)
530 profile1
[i
][j
] = (surf
[0][n
][j
+ ybins
* i
]);
531 profile2
[i
][j
] = (surf
[1][n
][j
+ ybins
* i
]);
532 /*Finding max and min values*/
533 if (profile1
[i
][j
] > max1
)
535 max1
= profile1
[i
][j
];
537 if (profile1
[i
][j
] < min1
)
539 min1
= profile1
[i
][j
];
541 if (profile2
[i
][j
] > max2
)
543 max2
= profile2
[i
][j
];
545 if (profile2
[i
][j
] < min2
)
547 min2
= profile2
[i
][j
];
552 write_xpm(xpmfile1
, 3, numbuf
, "Height", "x[nm]", "y[nm]", xbins
, ybins
, xticks
, yticks
,
553 profile1
, min1
, max1
, lo
, hi
, &maplevels
);
554 write_xpm(xpmfile2
, 3, numbuf
, "Height", "x[nm]", "y[nm]", xbins
, ybins
, xticks
, yticks
,
555 profile2
, min2
, max2
, lo
, hi
, &maplevels
);
558 gmx_ffclose(xpmfile1
);
559 gmx_ffclose(xpmfile2
);
569 static void writeraw(real
*** surf
, int tblocks
, int xbins
, int ybins
, gmx::ArrayRef
<const std::string
> fnms
)
574 raw1
= gmx_ffopen(fnms
[0], "w");
575 raw2
= gmx_ffopen(fnms
[1], "w");
576 fprintf(raw1
, "#Legend\n#TBlock\n#Xbin Ybin Z t\n");
577 fprintf(raw2
, "#Legend\n#TBlock\n#Xbin Ybin Z t\n");
578 for (n
= 0; n
< tblocks
; n
++)
580 fprintf(raw1
, "%5d\n", n
);
581 fprintf(raw2
, "%5d\n", n
);
582 for (i
= 0; i
< xbins
; i
++)
584 for (j
= 0; j
< ybins
; j
++)
586 fprintf(raw1
, "%i %i %8.5f\n", i
, j
, (surf
[0][n
][j
+ ybins
* i
]));
587 fprintf(raw2
, "%i %i %8.5f\n", i
, j
, (surf
[1][n
][j
+ ybins
* i
]));
597 int gmx_hydorder(int argc
, char* argv
[])
599 static const char* desc
[] = {
600 "[THISMODULE] computes the tetrahedrality order parameters around a ",
601 "given atom. Both angle an distance order parameters are calculated. See",
602 "P.-L. Chau and A.J. Hardwick, Mol. Phys., 93, (1998), 511-518.",
603 "for more details.[PAR]",
604 "[THISMODULE] calculates the order parameter in a 3d-mesh in the box, and",
605 "with 2 phases in the box gives the user the option to define a 2D interface in time",
606 "separating the faces by specifying parameters [TT]-sgang1[tt] and",
607 "[TT]-sgang2[tt] (it is important to select these judiciously)."
611 static int nsttblock
= 1;
612 static int nlevels
= 100;
613 static real binwidth
= 1.0; /* binwidth in mesh */
615 static real sg2
= 1; /* order parameters for bulk phases */
616 static gmx_bool bFourier
= FALSE
;
617 static gmx_bool bRawOut
= FALSE
;
618 int frames
, xslices
, yslices
; /* Dimensions of interface arrays*/
619 real
*** intfpos
; /* Interface arrays (intfnr,t,xy) -potentially large */
620 static const char* normal_axis
[] = { nullptr, "z", "x", "y", nullptr };
623 { "-d", FALSE
, etENUM
, { normal_axis
}, "Direction of the normal on the membrane" },
624 { "-bw", FALSE
, etREAL
, { &binwidth
}, "Binwidth of box mesh" },
625 { "-sgang1", FALSE
, etREAL
, { &sg1
}, "tetrahedral angle parameter in Phase 1 (bulk)" },
626 { "-sgang2", FALSE
, etREAL
, { &sg2
}, "tetrahedral angle parameter in Phase 2 (bulk)" },
627 { "-tblock", FALSE
, etINT
, { &nsttblock
}, "Number of frames in one time-block average" },
628 { "-nlevel", FALSE
, etINT
, { &nlevels
}, "Number of Height levels in 2D - XPixMaps" }
632 /* files for g_order */
633 { efTRX
, "-f", nullptr, ffREAD
}, /* trajectory file */
634 { efNDX
, "-n", nullptr, ffREAD
}, /* index file */
635 { efTPR
, "-s", nullptr, ffREAD
}, /* topology file */
636 { efXPM
, "-o", "intf", ffWRMULT
}, /* XPM- surface maps */
637 { efOUT
, "-or", "raw", ffOPTWRMULT
}, /* xvgr output file */
638 { efOUT
, "-Spect", "intfspect", ffOPTWRMULT
}, /* Fourier spectrum interfaces */
640 #define NFILE asize(fnm)
643 const char * ndxfnm
, *tpsfnm
, *trxfnm
;
644 gmx_output_env_t
* oenv
;
646 if (!parse_common_args(&argc
, argv
, PCA_CAN_VIEW
| PCA_CAN_TIME
, NFILE
, fnm
, asize(pa
), pa
,
647 asize(desc
), desc
, 0, nullptr, &oenv
))
651 bFourier
= opt2bSet("-Spect", NFILE
, fnm
);
652 bRawOut
= opt2bSet("-or", NFILE
, fnm
);
656 gmx_fatal(FARGS
, "Can not have binwidth < 0");
659 ndxfnm
= ftp2fn(efNDX
, NFILE
, fnm
);
660 tpsfnm
= ftp2fn(efTPR
, NFILE
, fnm
);
661 trxfnm
= ftp2fn(efTRX
, NFILE
, fnm
);
664 GMX_RELEASE_ASSERT(normal_axis
[0] != nullptr,
665 "Option setting inconsistency; normal_axis[0] is NULL");
666 if (std::strcmp(normal_axis
[0], "x") == 0)
670 else if (std::strcmp(normal_axis
[0], "y") == 0)
674 else if (std::strcmp(normal_axis
[0], "z") == 0)
680 gmx_fatal(FARGS
, "Invalid axis, use x, y or z");
685 case 0: fprintf(stderr
, "Taking x axis as normal to the membrane\n"); break;
686 case 1: fprintf(stderr
, "Taking y axis as normal to the membrane\n"); break;
687 case 2: fprintf(stderr
, "Taking z axis as normal to the membrane\n"); break;
690 /* tetraheder order parameter */
691 /* If either of the options is set we compute both */
692 gmx::ArrayRef
<const std::string
> intfn
= opt2fns("-o", NFILE
, fnm
);
693 if (intfn
.size() != 2)
695 gmx_fatal(FARGS
, "No or not correct number (2) of output-files: %td", intfn
.ssize());
697 calc_tetra_order_interface(ndxfnm
, tpsfnm
, trxfnm
, binwidth
, nsttblock
, &frames
, &xslices
,
698 &yslices
, sg1
, sg2
, &intfpos
, oenv
);
699 writesurftoxpms(intfpos
, frames
, xslices
, yslices
, binwidth
, intfn
, nlevels
);
703 gmx::ArrayRef
<const std::string
> spectra
= opt2fns("-Spect", NFILE
, fnm
);
704 if (spectra
.size() != 2)
706 gmx_fatal(FARGS
, "No or not correct number (2) of output-files: %td", spectra
.ssize());
708 powerspectavg(intfpos
, frames
, xslices
, yslices
, spectra
);
713 gmx::ArrayRef
<const std::string
> raw
= opt2fns("-or", NFILE
, fnm
);
716 gmx_fatal(FARGS
, "No or not correct number (2) of output-files: %td", raw
.ssize());
718 writeraw(intfpos
, frames
, xslices
, yslices
, raw
);