The crude method would be to sort() the data, and divide up the resulting vector into equal chunks and look at the values there to determine appropriate boundaries.
This obviously has difficulties with duplicate entries that happen to fall at the boundaries. I have no idea how that is normally dealt with. Consider for example that bimodal distributions (only two distinct values) cannot be divided up into distinct quartiles, and even if you toss in a small number of other values so that you have at least as many distinct values as quantiles, clearly you cannot promise that you will get bins of the same width.
But sort and divide up might be good enough for your purposes.
