data {
   for(i in 1:N) {
      for(j in 1:T) {
         ## risk set = 1 if obs.t >= t
         Y[i,j] <- step(obs.t[i] - t[j] + eps)
                                                                                
         ## counting process jump = 1 if obs.t in [ t[j], t[j+1] )
         ##                      i.e. if t[j] <= obs.t < t[j+1]
         dN[i,j] <- Y[i,j] * fail[i] * step(t[j+1] - obs.t[i] - eps)
      }
   }
}
model {
   for(j in 1:T) {
      for(i in 1:N) {
         dN[i,j]   ~ dpois(Idt[i,j]);              # Likelihood
         Idt[i,j] <- Y[i,j] * exp(beta*Z[i]) * dL0[j]; # Intensity 
      }     
      dL0[j] ~ dgamma(mu[j], c);

      mu[j] <- dL0.star[j] * c;    # prior mean hazard
      dL0.star[j] <- r * (t[j+1]-t[j])  
   
      # Survivor function = exp(-Integral{l0(u)du})^exp(beta*z)    
      S.treat[j] <- pow(exp(-sum(dL0[1:j])), exp(beta * -0.5));
      S.placebo[j] <- pow(exp(-sum(dL0[1:j])), exp(beta * 0.5));	
   }
   beta ~ dnorm(0.0, 1.0E-6);                 
    
   c <- 0.001; 
   r <- 0.1; 
}