data {

  int N;  // # observations
  int J;  // # chicks

  real weight[N];
  real Time[N];

  int Chick[N]; // needs to be integers, not a factor variable

}
  
transformed data {

  real y[N];

  y <- log(weight);

}
  
parameters{

  real b0;

  real b1;

  real eta0[J];

  real eta1[J];

  real<lower=0> sigma;

  real<lower=0> tau0;

  real<lower=0> tau1;

}
  
transformed parameters {

  real mu[N];

  real a0[J];

  real a1[J];

  for (j in 1:J) {

    a0[j] <- b0 + eta0[j];

    a1[j] <- b1 + eta1[j];

  }

  for (i in 1:N) {

    mu[i] <- a0[Chick[i]] + a1[Chick[i]]*Time[i];
  
  }

}
  
model {

  // Level 1

  for (i in 1:N) {

    y[i] ~ normal(mu[i],sigma) ;

  }

  // Level 2

  for (j in 1:J) {

    eta0[j] ~ normal(0,tau0);

    eta1[j] ~ normal(0,tau1);

  }

  // Prior distributions

  b0 ~ normal(0,1e6);

  b1 ~ normal(0,1e6);

  sigma ~ uniform(0,50);

  tau0 ~ uniform(0,50);

  tau1 ~ uniform(0,50);

}
  
generated quantities {

  real yrep[N];

  for (i in 1:N) {

    yrep[i] <- normal_rng(mu[i],sigma);

  }

}

CXXFLAGS = -O2 -Wall $(DEBUGFLAG) -mfpmath=sse -msse2 -mstackrealign $(LTO)