I reran the numbers - disclaimer, it has been a long while since I've done this... easier to run in Matlab if someone was so inclined...
Let's say we have a 16 bit dac with LSB=100 fF (σ=0.1%) and an MSB = 2^15*100fF = 3.2768 nF (σ=5.5e-6). The σ of the MSB (8000) will be 18 fF. Also, the σ of the 7FFF code will be essentially 18 fF since the cap for that code is only 100 fF less. Thus the σ of the difference will be sqrt(2)*18 fF = 25 fF, or a quarter LSB. This will give you the DNL σ so if you are only interested in DNL (or monotonicity) it is better than 16 bits (one σ). I think this is how your formula analyzes it (but perhaps it should be 2
N-1*1u for looking at the quantization noise at the major transition).
In hindsight that isn't so surprising since 100 fF is a rather large LSB cap resulting in a total cap of 6.4 nF.
But INL is what causes the big tones/intermodulation. INL σ isn't straight forward - I forgot the formula (which is only an approximation anyway) but I think the worst case is on the order of sqrt(2^16)*σ
LSB which is 256 LSBs. INL will give you tones depending on the shape so figuring the tone σ isn't straightforward, but it is going to be closer to 8-10 bits than 16 bits. However, I believe that formula overestimated the INL so you really have to run a Monte Carlo or find someone with a better memory than me
.