It seems that collect(f,[xa^2,xb^2,xc^2]) does not work as you would desire, or as I would expect.
A work around is the command:
This substitutes instances of squares to new variables, collects in these variables, then substitutes back to the squares of the original variables.
This is not elegant, and will not work when your polynomial involves third powers or more. It will work for your particular example. You could generalize this to higher degrees by using algsubs (which will convert xa^3 to a*xa etc.), though you'd need to make a version of algsubs which can take a list of substitutions.
Note that collect here first collects it as a polynomial in xa^2, then writes the coefficients as polynomials in xb^2, etc., recursively.