This one weird trick integrates fractals. But does it deliver the correct results?

There is a wonderful paper in the April 2000 issue of the Monthly by Bob Strichartz entitled "Evaluating Integrals Using Self-Similarity". The MAA has PDF copy here. You've effectively got the basics of this technique, though, there is much more that one can do.

In general, any self-similar set supports a multitude of self-similar measures and you can use the self-similarity to integrate with respect to those measures. A recursive procedure, like yours, can be used to compute the integral of any polynomial. We can also use the self-similarity to generate a numerical procedure to estimate integrals of more arbitrary functions.

As you note, the unit interval is self-similar and it's not too hard use this to prove that $$\int_0^1 x^n \, dx = \frac{1}{n+1}.$$ However, other ways of decomposing the interval into intervals of different sizes or just weighting the subparts different amounts leads to other self-similar measures on the interval.

And, yes, your computation is certainly correct - as is Jack's.

It's all really quite beautiful stuff.

It looks right to me. You may notice that there is an even faster way to evaluate the integral: since the Canton set is symmetric with respect to $x=\frac{1}{2}$, $$ I=\int x\,d\mu_C = \frac{1}{2}\left(\int x\,d\mu_C+\int (1-x)\,d\mu_C\right)=\frac{1}{2}\int 1\,d\mu_C = \frac{1}{2}.$$ The Hausdorff dimension does not matter, we are simply putting a uniform measure over a measurable subset of $\mathbb{R}$. Moreover, by merging this technique and yours, it is not difficult to check that: $$ \int x^2\,d\mu_C = \frac{3}{8}.$$