With the question by the OP asking for a reference I will try to do
just that, providing links to the OEIS. Surprisingly enough even the
OEIS does not offer the usual variety of references here, suggesting
that this problem is open. We compute generating functions $T_{\le
h}(z)$ for the height being at most $h$ and the desired count is then
given by $T_{\le h}(z)-T_{\le h-1}(z).$
We have for the labeled the combinatorial class
$$\mathcal{T}_{\le h} =
\mathcal{Z}\def\textsc#1{\dosc#1\csod}
\def\dosc#1#2\csod{{\rm #1{\small #2}}}\times \textsc{SET}(\mathcal{T}_{\le h-1}).$$
For the labeled case this translates into
$$T_{\le h}(z) = z\exp T_{\le h-1}(z)$$
with $T_{\le 0}(z) = 0$ and $T_{\le 1}(z) = z.$
For example with $h=3$ we get
$$T_{\le 3}(z) = z \exp(z \exp(z))$$
which yields the sequence OEIS A052512:
$$1, 2, 9, 40, 205, 1176, 7399, 50576, \ldots$$
and for $h=4$ we get
$$T_{\le 4}(z) = z(\exp( z \exp(z \exp(z)))$$
which yields the sequence OEIS A052513:
$$1, 2, 9, 64, 505, 4536, 46249, 526352, \ldots$$
We have for the unlabeled the combinatorial class
$$\mathcal{T}_{\le h} =
\mathcal{Z}\times \textsc{MSET}(\mathcal{T}_{\le h-1}).$$
which translates into generating functions, as follows,
$$T_{\le h}(z) = z\exp
\left(\sum_{l\ge 1}
\frac{T_{\le h-1}(z^l)}{l}\right).$$
where $T_{\le 0}(z) = 0$ and $T_{\le 1}(z) = z$ as before.
We get for
$$T_{\le 2}(z) =
z\exp \left(\sum_{l\ge 1} \frac{z^l}{l}\right)
= z\exp\log\frac{1}{1-z} = \frac{z}{1-z}.$$
These are root nodes with singletons attached to them.
Continuing with the recursion we obtain
$$T_{\le 3}(z) = z\exp
\left(\sum_{l\ge 1}
\frac{1}{l}\sum_{q\ge 1} z^{ql}\right)
= z\exp
\left(\sum_{q\ge 1}
\sum_{l\ge 1} \frac{1}{l} z^{ql}\right)
\\ = z\exp\sum_{q\ge 1}\log\frac{1}{1-z^q}
= z \prod_{q\ge 1} \frac{1}{1-z^q}.$$
which yields the sequence OEIS A000041:
$$1, 1, 2, 3, 5, 7, 11, 15, 22, 30, 42, 56, \ldots$$
(partition numbers, not surprisingly, as we allocate one node for the
root and the rest is a partition into singleton fans of height at most
two, of which there is only one per node count).
From $T_{\le 4}(z)$ we obtain the sequence
OEIS A001383:
$$1, 1, 2, 4, 8, 15, 29, 53, 98, 177, 319, 565, \ldots$$
Finally $T_{\le 5}(z)$ yields the sequence
OEIS A001384:
$$1, 1, 2, 4, 9, 19, 42, 89, 191, 402, 847 \ldots$$
There are two interpretations here depending on whether the singleton
is supposed to have height zero or one. A relevant link is this
MSE link.
Remark. The complexity of the growth of the partition numbers
provides an idea of the difficulty of this problem, which is as
difficult if not more. Studying the OEIS entries in more detail it
does appear that simple recurrences for these can be computed.
Concerning the recurrence for the labeled case differentiation
produces
$$T'_{\le h}(z) = \exp T_{\le h-1}(z) +
z \exp T_{\le h-1}(z) \times T'_{\le h-1}(z)
\\ = \frac{1}{z} T_{\le h}(z) + T_{\le h}(z) \times T'_{\le h-1}(z).$$
Extracting coefficients here we get
$$n! [z^n] T'_{\le h}(z)
= T_{\le h, n+1} = n! [z^{n+1}] T_{\le h}(z)
+ n! [z^n ] T_{\le h}(z) \times T'_{\le h-1}(z)
\\ = \frac{1}{n+1} T_{\le h, n+1}
+ n! \sum_{q=1}^{n} [z^q] T_{\le h}(z) [z^{n-q}] T'_{\le h-1}(z)
\\ = \frac{1}{n+1} T_{\le h, n+1}
+ \sum_{q=1}^{n} {n\choose q} T_{\le h,q} T_{\le h-1, n-q+1}.$$
We obtain the closed form
$$\bbox[5px,border:2px solid #00A000]
{T_{\le h, n+1} = \frac{n+1}{n}
\sum_{q=1}^{n} {n\choose q} T_{\le h,q} T_{\le h-1, n-q+1}.}$$
The boundary condition here is $T_{\le h, 1} = 1$ and
$T_{0,1} = 0.$
Applying differentiation to the unlabeled case yields
$$T'_{\le h}(z) =
\exp\left(\sum_{l\ge 1} \frac{T_{\le h-1}(z^l)}{l}\right)
\\ + z\exp \left(\sum_{l\ge 1} \frac{T_{\le h-1}(z^l)}{l}\right)
\times \sum_{l\ge 1} T'_{\le h-1}(z^l) z^{l-1}
\\ = \frac{1}{z} T_{\le h}(z)
+ T_{\le h}(z) \times \sum_{l\ge 1} T'_{\le h-1}(z^l) z^{l-1}.$$
Extracting coefficients we find
$$[z^n] T'_{\le h}(z) = (n+1) T_{\le h, n+1}
\\ = [z^{n+1}] T_{\le h}(z)
+ \sum_{l=1}^n
\sum_{q=1}^n
[z^q] T_{\le h}(z) [z^{n-q}] z^{l-1} T'_{\le h-1}(z^l)
\\ = T_{\le h, n+1}
+ \sum_{l=1}^n
\sum_{q=1}^n
T_{\le h,q} [z^{n-q-l+1}] T'_{\le h-1}(z^l).$$
We obtain the closed form
$$T_{\le h, n+1} =
\frac{1}{n} \sum_{l=1}^n
\sum_{q=1}^n
T_{\le h,q} [[l|n-q+1]] ((n-q+1)/l) T_{\le h-1, (n-q+1)/l}$$
which simplifies to
$$\bbox[5px,border:2px solid #00A000]
{T_{\le h, n+1} =
\frac{1}{n}
\sum_{q=1}^n \sum_{l|n-q+1}
T_{\le h,q} ((n-q+1)/l) T_{\le h-1, (n-q+1)/l}.}$$
Boundary conditions are the same as in the labeled case.
This recurrence makes it possible to compute e.g. $T_{\le 6,n}$
which yields
$$1, 1, 2, 4, 9, 20, 47, 108, 252, 582,
\\ 1345, 3086, 7072, 16121,\ldots$$
which is OEIS A001385 and $T_{\le 7,n}$
which yields
$$1, 1, 2, 4, 9, 20, 48, 114, 278, 676,
\\ 1653, 4027, 9816, 23843,\ldots$$
which is OEIS A034823.
The Maple code to implement these two recurrences is as follows.
T :=
proc(h,m)
option remember;
local n;
if m=1 then if h=0 then return 0 else return 1 fi fi;
n := m-1;
(n+1)/n*
add(binomial(n,q)*T(h,q)*T(h-1,n-q+1), q=1..n);
end;
X :=
proc(h,m)
option remember;
local n;
if m=1 then if h=0 then return 0 else return 1 fi fi;
n := m-1;
1/n*
add(add(X(h,q)*
`if`(n-q+1 mod l =0, (n-q+1)/l*X(h-1,(n-q+1)/l), 0),
q=1..n), l=1..n);
end;
Best Answer
The number of rooted, ordered, incomplete, unlabeled $k$-ary trees with $n$ vertices is given by
$$C^{(k)}_n=\frac1{(k-1)n+1}\binom{kn}n\;.$$
These are sometimes called Fuss-Catalan numbers; see Concrete Mathematics (p. 347) and MathWorld (which gives two references). Their generating function $C^{(k)}(x)=\sum_0^\infty C^{(k)}_nx^n$ satisfies $C^{(k)}(x)=1+xC^{(k)}(x)^k$. The numbers of rooted, ordered, incomplete, unlabeled ternary ($k=3$), quartic ($k=4$), qunitic ($k=5$), sextic ($k=6$), heptic ($k=7$) and octic ($k=8$) trees form OEIS sequences A001764, A002293, A002294, A002295, A002296 and A007556, respectively. To get the number of labeled trees, just multiply by $n!$.