*
* @ARDLM(options) LagVect
*
* Generates the state-space system matrices for a univariate
* autoregressive process whose lag coefficients are in the parameter
* LagVect VECTOR. You can also use it to set up the framework (proper
* dimensions, all the fixed values put in place) based upon an input
* number of lags (LAGS option) to fill out the non-fixed elements later.
*
* The simplest state-space (dynamic linear model) representation of
*
*  y(t)=a1 y(t-1) + a2 y(t-2) + ... + ap y(t-p) + u(t)
*
* has state-vector X(t)=[y(t),y(t-1),...,y(t-p+1] and
*
*     |a1,a2,...,ap-1,ap|
*     |1 , 0,...,0   ,0 |
* A = |0 , 1,...,0   ,0 |
*     | ...     ,0   ,0 |
*     |0 , 0,...,1   ,0 |
*
* F = p dimension unit vector
*
* C = p dimension unit vector
*
* Note that the only thing that depends upon the free coefficients is
* the first row in A. In this form, there is no "measurement" error as
* the only shock (u) enters the state equation.
*
* Parameters:
*
*  LagVect (optional, input)  VECTOR of AR coefficients on consecutive lags
*
*  Options:
*
*  LAGS=number of lags in model (input)
*    Can be used as an alternative to the LagVect parameter to set up the
*    framework of the matrices if you don't yet know the coefficients.
*
*  A=DLM system "A" matrix (output; input/output if used with POKE)
*  F=DLM system "F" matrix (output)
*  C=DLM system "C" matrix (output)
*
*  POKE/[NOPOKE]
*    If POKE, assumes that the A option matrix has the framework set up and
*    only changes the first row to insert the coefficients out of LagVect.
*    This avoids repeating the set up of the fixed elements.
*
* Revision Schedule:
*  07/2025 Written by Tom Doan, Estima
*
procedure ARDLM LagVect
type vect LagVect
*
option rect *A
option rect *F
option rect *C
option int  lags
option switch poke 0
*
local int workLags stateLags inputLags l
*
* Figure out for how many lags to set up the model. This is either the
* dimension of LagVect or the value of the LAGS option.
*
if %defined(LagVect)
   compute inputLags=%size(LagVect)
else
   compute inputLags=0
*
compute workLags=%imax(inputLags,lags)
*
* This is the size of the state vector, which has to be at least one.
*
compute stateLags=%imax(workLags,1)
*
compute dim=stateLags
if %defined(A) {
   if .not.poke
      compute A=%zeros(dim,dim)
   do l=1,workLags
      if l<=inputLags
         compute A(1,l)=LagVect(l)
      if l<>1
         compute A(l,l-1)=1
   end do l
}
if %defined(F) {
   compute F=%zeros(dim,1)
   compute F(1,1)=1.0
}
if %defined(C) {
   compute C=%zeros(dim,1)
   compute C(1,1)=1.0
}
end
