Taking advantage of the spherical symmetry to decompose this into a sum of one-dimensional problems sounds like the right approach. I will probably just be redoing what Reed-Simon had in mind here.

The Laplacian
is given in polar coordinates by
$$
\Delta=\frac{1}{r}\frac{\partial}{\partial r} r\frac{\partial}{\partial r} + \frac{1}{r^2}\frac{\partial^2}{\partial\theta^2} .
$$
We identify $L^2(\mathbb R^2)=L^2((0,\infty); rdr)\otimes L^2(S^1)$ and decompose $L^2(S^1)=\bigoplus L(e^{in\theta})$ and also get rid of the weight in the radial component with the help of the unitary map $U:L^2(r\, dr)\to L^2(0,\infty)$, $(Uf)(r)=r^{1/2}f(r)$, then $H$ becomes (is unitarily equivalent to, to be precise) an orthogonal sum $H=\bigoplus H_n$, with each $H_n$ acting in $L^2(0,\infty)$ and
$$
H_n f = -f''+\frac{n^2-1/4+c}{r^2}f .
$$
This operator is essentially self-adjoint if and only if $n^2-1/4+c\ge 3/4$ or $c\ge 1-n^2$. See for example the corresponding discussion in my answer here.

So $H$ on $C_0^{\infty}(\mathbb R^2\setminus\{ 0\})$ itself is essentially self-adjoint if and only if this happens for all $n\in\mathbb Z$, which is indeed equivalent to $c\ge 1$. For $0<c<1$, $H_n$ is still essentially self-adjoint on $C_0^{\infty}(0,\infty)$ for $|n|\ge 1$, while the closure of $H_0$ on this domain is symmetric with deficiency index $(1,1)$. So in this case there is a one-parameter family of boundary conditions, to be imposed on the elements of $D(H_0^*)$, and this gives all self-adjoint realizations.

As for (1), we can just set $c=0$ in this analysis and everything is still correct, and we are now discussing the self-adjoint extensions of $-\Delta$ on $C_0^{\infty}(\mathbb R^2\setminus \{ 0\})$. Since there is nothing to choose in the higher harmonics, the question is if $-d^2/dr^2+k/r^2$ can ever be given the same domain as $-d^2/dr^2-1/(4r^2)$ for $-1/4<k<3/4$. The discussion here strongly suggests that the answer is *no*, though it's not completely conclusive since we fixed specific boundary conditions there.

**Added later (with one detail corrected later still):** This seems quite clear actually: When $c=0$, the possible domains $D(H_0)$ can be obtained by putting an arbitrary solution of $-f''-f/(4r^2)=0$ into $D(H_0)$. This is an Euler equation with solution basis $f_1=r^{1/2}$, $f_2=r^{1/2}\log r$. Note also that when we undo the unitary transformation from above, these correspond to $U^*f_j=1, \log r$. On the other hand, $f$ can never be in $D(H_0)$ when $c>0$ because then $-f''+(c-1/4)f/r^2= cf/r^2\notin L^2$. So we never have $D(H)=D(\Delta)$.

For the record ("what is the correction?"), in an earlier version I put $f_1=r^{1/2}$ into $D(H_0)$ and claimed that this corresponds to $D_0=H_0^1\cap H^2$, which you suggested as a possible domain. This was nonsense of course (I mixed up my signs in the exponent when applying $U^*$), since in that case we put $g=1$ (near $r=0$) into $D(\Delta)$ and we simply obtain the plain Laplacian on $H^2(\mathbb R^2)$. Any other choice, namely put $g=b+\log r$ into $D(\Delta)$, is certainly closer to what you had in mind, but it's not really the domain $D_0$. In fact, $H_0^1$ on $\mathbb R^2\setminus \{0\}$ is simply $H^1$ since the function $h(r)=\log a^2/\log r^2$ gets us from $h(0)=0$ to $h(a)=1$ and $\|h\|_{H^1(r<a)}\to 0$ as $a\to 0$.

Finally, (2) is now moot, but it would be a consequence of (1) since if $D(S)\supseteq D(T)$ for closed operators, then $S$ is $T$-bounded (= Theorem 5.9 in Weidmann, Linear operators in Hilbert spaces).