Implementing bounded parameters in JAGS

Question

I have a simple hierarchical model that I am implementing in JAGS where I model length as a function of age as a single change-point regression, then model these variables against a site level parameter called temp. I would like to keep the cp parameter bounded within the minimum and maximum age for each site. I've tried looking at both dinterval and the T(,) construct in JAGS, but neither seem to do quite what I want. In my reprex below the dinterval function simply limits cp to 0 instead of the minimum age at a site.

EDIT: I realized after initial post that the parameters must be bounded on theta_hat, not on individual parameter naming as I originally posted (i.e. cp[j] <- theta_hat[j, 4]). Otherwise, the second-level model predictions will be incorrect.

Any thoughts or ideas for how I can keep theta_hat parameters bounded?

library(R2jags)

age1 <- sort(rep(2:8, 10))
age2 <- sort(rep(1:10, 10))
vb <- function(age, linf, k, t0) linf * (1 - exp(-k * (age - t0)))
len1 <- vb(age1, 100, 0.25, -0.5) + rnorm(length(age1), 0, age1)
len2 <- vb(age2, 80, 0.5, -0.75) + rnorm(length(age2), 0, age2)
laa_data1 <- data.frame(site = 1L, age = age1, len = len1)
laa_data2 <- data.frame(site = 2L, age = age2, len = len2)
laa_data <- rbind(laa_data1, laa_data2)

K <- 4 # number of parameters
W <- diag(4)

## Model ##
data = list(
  dat = as.matrix(laa_data),
  n = nrow(laa_data),
  nsites = 2,
  temp = as.vector(c(20,25)),
  # minage = c(min(age1), min(age2)),
  # maxage = c(max(age1), max(age2)),
  K = K,
  W = W
)

model <- function() {

for (j in 1:nsites){
  
  alpha[j] <- theta[j, 1]
  beta1[j] <- theta[j, 2]
  beta2[j] <- theta[j, 3]
  cp[j] <-    theta[j, 4]
  

  theta[j, 1:K] ~ dmnorm(theta_hat[j, ], theta_tau[,]) # is there a way to use T(,) to keep parameters bounded (or censored) here?

  # or should parameter bounding happen at this level?
  theta_hat[j, 1] <- a1 + (b1 * temp[j])
  theta_hat[j, 2] <- a2 + (b2 * temp[j])
  theta_hat[j, 3] <- a3 + (b3 * temp[j])
  theta_hat[j, 4] <- a4 + (b4 * temp[j])
}

# likelihood
for (i in 1:n) {
  mu[i] <- 
    alpha[dat[i,1]] + 
    beta1[dat[i,1]] * min(dat[i,2], cp[dat[i,1]]) + 
    beta2[dat[i,1]] * max(0,dat[i,2] - cp[dat[i,1]])
  dat[i,3] ~ dnorm(mu[i], 1/nu_alk)
}

# priors
a1 ~ dunif(-10, 10)
b1 ~ dunif(-10, 10) 
a2 ~ dunif(-10, 10)
b2 ~ dunif(-10, 10)
a3 ~ dunif(-10, 10)
b3 ~ dunif(-10, 10)
a4 ~ dunif(-10, 10)
b4 ~ dunif(-10, 10)

# priors on variance-covariance
theta_tau[1:K,1:K] ~ dwish(W[,], df)
df <- K + 1
nu_alk ~ dunif(0, 100)

}

params = c("theta_hat", "alpha", "beta1", "beta2", "cp") 

# Fit model
fit <- jags(
  data=data,
  parameters.to.save = params,
  model=model, 
  n.chains = 1,
  n.iter = 5000, 
  n.burnin = 1000,
  n.thin=5
)