TL;DR: Yes "th" is gradually increasing and possibly exceeding bounds of (-pi,pi], but it doesn't matter. You'll always get the "right" quaternion.
What would cause the quaternion value be NaN? This is actually one of the benefits of a quaternion -- for describing any arbitrary rotation, there are only two possible quaternions. This is opposed to a 3-parameter description of rotation (like any of the family of Euler angles) where there are infinite choices that all yield the same rotation. In your specific question the same quaternion will be produced for $th = th + 2 n pi$ for any choice of integer n. In other words for any yaw angle you'll always get the same quaternion regardless of how many "wraps" around the unit circle you have.
If you look at the source code for createQuaternionMsgFromYaw, you'll see that they actually call tf::Quaternion.setRPY to calculate the values of the quaternion. If you look at the source code for the setRPY method you'll see it's just a bunch of sin and cos calls that also don't matter how many wraps you have.
