Problem with the derivative of Heaviside distribution
Naively we have $\frac{d}{dx}\theta(y(x)) = \theta'(y(x)) \, y'(x)$ by the chain rule.
By implicit derivation and chain rule we have $$ y'(x) = -\operatorname{sign}(x-f(y)) \, (1-f'(y)y'(x)), $$ which when solved for $y'(x)$ gives $$ y'(x) = \frac{\operatorname{sign}(x-f(y))}{\operatorname{sign}(x-f(y)) f'(y) - 1}. $$
Thus, since $\theta'=\delta,$ we get $$ \frac{d}{dx}\theta(y(x)) = \delta(y(x)) \, \frac{\operatorname{sign}(x-f(y))}{\operatorname{sign}(x-f(y)) f'(y) - 1} = \delta(y(x)) \, \frac{\operatorname{sign}(x-f(0))}{\operatorname{sign}(x-f(0)) f'(0) - 1}, $$ where it in the last step has been used that $g(x)\delta(x) = g(0)\delta(x).$
As I wrote at the start of this post, these calculations were done naively. But still I think that they are valid.
After looking a bit more into this I realized that since $$ \delta(y(x)) = \sum_{z \in y^{-1}(\{0\})} \frac{\delta(x-z)}{|y'(z)|} $$ we have $$ \frac{d}{dx}\theta(y(x)) = \theta(y(x)) \, y'(x) = \delta(y(x)) \, y'(x) \\ = \sum_{z \in y^{-1}(\{0\})} \frac{y'(x)}{|y'(z)|} \delta(x-z) \\ = \sum_{z \in y^{-1}(\{0\})} \operatorname{sign}(y'(z)) \delta(x-z) \\ = \sum_{z \in y^{-1}(\{0\})} \operatorname{sign}\left(\frac{\operatorname{sign}(z-f(y(z)))}{\operatorname{sign}(z-f(y(z))) f'(y(z)) - 1}\right) \delta(x-z) \\ = \sum_{z \in y^{-1}(\{0\})} \operatorname{sign}\left(\frac{\operatorname{sign}(z-f(0))}{\operatorname{sign}(z-f(0)) f'(0) - 1}\right) \delta(x-z) . $$