Hidden Markov Model: emissions probabilities dependent on observable parameter

I need to fit an HMM where the emission probabilities (ep) are discrete and dependent on a known variable quantity. E.g.:

Imagine a daily time series of binary emissions ("1" or "2"). I suppose an hidden markov chain of 3 states $s_i$, each with probability $p_{s_i}$ of emitting "1" and $1-p_{s_i}$ of emitting "2".

When $p$ depends only on the state $s$ or on other fixed parameters (i.e. if $p_s$ were to be a NormalDistribution with unknown but time-constant $\mu_s$ and $\sigma_s$ for instance) it is straightforward to use EstimatedProcess to find a transition matrix and the emission probabilities from a dataset.

But how to include a dependence on a parameter that depends on the actual day $d$ of the time series (thus something known)?

The only way I could think of was to fictitiously create bivariate emissions, of the form ("1", $d$) or ("2", $d$), and then multiply the number of states by the cardinality of $d$ (in my case I need 3 hidden states and 7 weekdays, so 21 states!), where only 3 states can emit each value of $d$ with probability 1. But then:

1) I would have to define a precise transition matrix, where each day-state can only go in the 3 nextday-states.

2) I would need the optimization process not to touch some of such transition probabilities and the fixed emission probabilities (i.e. $s_{i,d^*}$ with $i=1,2,3$ must emit only ("1", $d^*$) or ("2", $d^*$) )

3) I would find myself with a very complex 21 states chain, that is probably overcomplicated for my problem...

Any ideas?