I would answer this first from the mathematical perspective. Mathematically, there is a difference between a quantity with a unit, like
$$1\ \mathrm{m}$$
and just
$$1$$
. This notion was - to a primeval extent - understood all the way back to Euclid of Alexandria (ca. 300 BCE), the same one for whom Euclidean geometry is eponymous. He considered the notion of magnitude, which had different "kinds" - lengths, areas, and so forth, against ratio, which was between two magnitudes of the same kind. Today, we might say magnitudes are the first type of quantity, and ratios as "pure numbers" of the second type of quantity.
The crucial difference is that magnitudes don't function like ratios do: you can add two magnitudes of the same kind, e.g. $1\ \mathrm{m} + 5\ \mathrm{m} = 6\ \mathrm{m}$, but you cannot add two magnitudes of different kinds, e.g. $1\ \mathrm{m} + 5\ \mathrm{m}^2$ is nonsense. You can multiply two magnitudes of different kinds, or the same kind, though - but the magnitude will be of another kind, e.g. two lengths multiply to an area, which we nowadays would use to describe something like
$$(a\ \mathrm{m})(b\ \mathrm{m}) = ab\ \mathrm{m^2}$$
But on the other hand, you can indeed add and multiply ratios as much as you want, because they're all the same type.
From the modern perspective, we would define ratios to be elements of the "pure" real number system $\mathbb{R}$, and while Euclid took magnitudes as his primitive notion from which are derived ratios, we take effectively ratios to be primary, and derive magnitudes from them. In a sense, elements of $\mathbb{R}$ are things that "want" to change the amount of something, to make it larger or smaller, or more or less numerous, according to the "sizing amount" they embody - if you follow this far enough you get to endomorphism rings and, even further, to the categorical product notion (see, e.g. Qiaochu Yuan's discussion here: https://math.stackexchange.com/questions/56663/is-there-a-natural-way-to-extend-repeated-exponentiation-beyond-integers/56710#56710). In fact, this is embodied in mathematical terminology: it's why we call them "scalars" in vector algebra. And it also connects with how we use them in natural language - when I say "three apples", "three" is something which "modifies", or acts upon, "apple" to conceptually produce the notion of three of them. Moreover, this is what multiplication is "really" about (instead of "repeated addition"): multiplication is, given two scalars, finding a third which acts upon things as one followed by another.
Magnitudes are a more complicated beast - the way you'd build them would be have to be basically a set of sets of "tagged" or "typed" real numbers that are effectively 1D vector spaces in that within any given such, tagged numbers can be added together and also rescaled by "pure" scalar (i.e. ratios) reals with an additional division operation which produces a ratio, and then finally a "trans-typal" multiplication operation that operates on the disjoint union of all the sets in the family (i.e. the "pooling" of all tagged reals). I am not sure at all what this kind of structure such a thing is called or if it's been studied in the mathematical literature before.
Now, returning to your question of physics and special relativity, we can employ this formalism to resolve the conundrum as thus. Instead of concentrating on all equations, we'll just concentrate on the basic one that describes the Minkowskian geometry, i.e. its line element
$$d\tau^2 = dt^2 - \frac{1}{c^2} (dx^2 + dy^2 + dz^2)$$
. Because we are talking physical quantities, these are magnitudes. The $dx$, etc. are NOT real numbers(*). They belong in the second kind of structure I just discussed. They have "units", and the typal multiplication by $\frac{1}{c^2}$ ends up converting units of (squared) space to units of (squared) time, so that the subtraction from $dt^2$ can proceed.
If, however, we take "$c = 1$", then technically there is a difference: for one, as written, that would be nonsense in the above as that would make $c$ a scalar ("ratio" in Euclidean terminology). As you note, we have to take $c = 1\ \frac{\mbox{distance unit}}{\mbox{time unit}}$, in order to stay within the confines of the magnitude system. However, if we do this, we can then pass through the natural isomorphism which "strips off the units", that is, which "pulls the tags from the tagged reals" that make up magnitudes, which we may denote as $\mathrm{strip}(x)$, and come to
$$\mathrm{strip}(d\tau^2) = \mathrm{strip}\left(dt^2 - \frac{1}{c^2}(dx^2 + dy^2 + dz^2)\right)$$
and then percolating all the isomorphism relations through,
$$[\mathrm{strip}(d\tau)]^2 = [\mathrm{strip}(dt)]^2 - \frac{1}{[\mathrm{strip}(c)]^2}([\mathrm{strip}(dx)]^2 + [\mathrm{strip}(dy)]^2 + [\mathrm{strip}(dz)]^2)$$
and since $\mathrm{strip}(c)$ is now a real number $1$ (take "$\frac{\mbox{distance unit}}{\mbox{time unit}}$" off), $\frac{1}{[\mathrm{strip}(c)]^2} = 1$ and thus if we abuse the notation on the $dx$, etc. values to go to their stripped counterparts, we get
$$d\tau^2 = dt^2 - (dx^2 + dy^2 + dz^2)$$
and it is in this sense that setting $c = 1$ makes things "unitless". If we did not do this, we could still pass through the $\mathrm{strip}$ isomorphism, but we'd end up with a factor $\frac{1}{c^2}$ which depended on the unit system we originally were coming from.
(*) Yes I'm calling a differential 1-form a real number - blah! I have to keep it simple somehow! Just pretend they're small, finite real numbers instead and the form is merely suggestive :)