pomp
Inference for partially observed Markov processes
All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
probe.c File Reference
#include "internal.h"
Include dependency graph for probe.c:

Go to the source code of this file.

Functions

SEXP apply_probe_data (SEXP object, SEXP probes)
 
SEXP apply_probe_sim (SEXP object, SEXP nsim, SEXP params, SEXP probes, SEXP datval, SEXP gnsi)
 

Function Documentation

◆ apply_probe_data()

SEXP apply_probe_data ( SEXP  object,
SEXP  probes 
)

Definition at line 4 of file probe.c.

4  {
5  SEXP retval, data, vals;
6  int nprobe;
7  int i;
8 
9  nprobe = LENGTH(probes);
10  PROTECT(data = GET_SLOT(object,install("data")));
11  PROTECT(vals = NEW_LIST(nprobe));
12  SET_NAMES(vals,GET_NAMES(probes));
13 
14  for (i = 0; i < nprobe; i++) {
15  SET_ELEMENT(vals,i,eval(PROTECT(lang2(VECTOR_ELT(probes,i),data)),
16  CLOENV(VECTOR_ELT(probes,i))));
17  if (!IS_NUMERIC(VECTOR_ELT(vals,i))) {
18  err("probe %d returns a non-numeric result",i+1);
19  }
20  UNPROTECT(1);
21  }
22  PROTECT(vals = VectorToPairList(vals));
23  PROTECT(retval = eval(PROTECT(LCONS(install("c"),vals)),R_BaseEnv));
24 
25  UNPROTECT(5);
26  return retval;
27 }
#define err(...)
Definition: pomp.h:21
Here is the caller graph for this function:

◆ apply_probe_sim()

SEXP apply_probe_sim ( SEXP  object,
SEXP  nsim,
SEXP  params,
SEXP  probes,
SEXP  datval,
SEXP  gnsi 
)

Definition at line 29 of file probe.c.

30  {
31  SEXP x, y, names, sims;
32  SEXP returntype, retval, val, valnames;
33  int nprobe, nsims, nobs, ntimes, nvals;
34  int xdim[2];
35  double *xp, *yp;
36  int p, s, i, j, k, len0 = 0, len = 0;
37 
38  PROTECT(nsim = AS_INTEGER(nsim));
39  if ((LENGTH(nsim)!=1) || (INTEGER(nsim)[0]<=0))
40  err("'nsim' must be a positive integer."); // #nocov
41 
42  PROTECT(gnsi = duplicate(gnsi));
43 
44  // 'names' holds the names of the probe values
45  // we get these from a previous call to 'apply_probe_data'
46  nprobe = LENGTH(probes);
47  nvals = LENGTH(datval);
48  PROTECT(names = GET_NAMES(datval));
49  PROTECT(returntype = NEW_INTEGER(1));
50  *(INTEGER(returntype)) = 0;
51  PROTECT(sims = do_simulate(object,params,nsim,returntype,gnsi));
52  PROTECT(y = VECTOR_ELT(sims,1));
53  *(INTEGER(gnsi)) = 0;
54 
55  nobs = INTEGER(GET_DIM(y))[0];
56  nsims = INTEGER(GET_DIM(y))[1];
57  ntimes = INTEGER(GET_DIM(y))[2];
58  // set up temporary storage
59  xdim[0] = nobs; xdim[1] = ntimes;
60  PROTECT(x = makearray(2,xdim));
61  setrownames(x,GET_ROWNAMES(GET_DIMNAMES(y)),2);
62 
63  // set up matrix to hold results
64  xdim[0] = nsims; xdim[1] = nvals;
65  PROTECT(retval = makearray(2,xdim));
66  PROTECT(valnames = NEW_LIST(2));
67  SET_ELEMENT(valnames,1,names); // set column names
68  SET_DIMNAMES(retval,valnames);
69 
70  for (p = 0, k = 0; p < nprobe; p++, k += len) { // loop over probes
71 
72  R_CheckUserInterrupt();
73 
74  for (s = 0; s < nsims; s++) { // loop over simulations
75 
76  // copy the data from y[,s,] to x[,]
77  xp = REAL(x);
78  yp = REAL(y)+nobs*s;
79  for (j = 0; j < ntimes; j++, yp += nobs*nsims) {
80  for (i = 0; i < nobs; i++, xp++) *xp = yp[i];
81  }
82 
83  // evaluate the probe on the simulated data
84  PROTECT(val = eval(PROTECT(lang2(VECTOR_ELT(probes,p),x)),
85  CLOENV(VECTOR_ELT(probes,p))));
86  if (!IS_NUMERIC(val)) {
87  err("probe %d returns a non-numeric result.",p+1);
88  }
89 
90  len = LENGTH(val);
91  if (s == 0)
92  len0 = len;
93  else if (len != len0) {
94  err("variable-sized results returned by probe %d.",p+1);
95  }
96  if (k+len > nvals)
97  err("probes return different number of values on different datasets.");
98 
99  xp = REAL(retval); yp = REAL(val);
100  for (i = 0; i < len; i++) xp[s+nsims*(i+k)] = yp[i];
101 
102  UNPROTECT(2);
103  }
104 
105  }
106  if (k != nvals)
107  err("probes return different number of values on different datasets."); // #nocov
108 
109  UNPROTECT(9);
110  return retval;
111 }
SEXP do_simulate(SEXP, SEXP, SEXP, SEXP, SEXP)
Definition: simulate.c:6
static R_INLINE void setrownames(SEXP x, SEXP names, int rank)
Definition: pomp_defines.h:110
static R_INLINE SEXP makearray(int rank, const int *dim)
Definition: pomp_defines.h:40
SEXP params
Definition: trajectory.c:128
Here is the call graph for this function:
Here is the caller graph for this function: