We present the consistent treatment of heavy quarks in jet production processes in DIS at NLO accuracy. The method is based on the ACOT massive factorization scheme (Aivazis-Collins-Olness-Tung) which is already in use in global fits for PDFs. Extension of this method to jets faces the problem of infra-red singularities, which have to cancel between real and virtual corrections. Moreover, potential mass singularities have to be controlled in both contributions. To this end, we extended the dipole subtraction method (S. Catani, S. Dittmaier, M. H. Seymour, Z. Trocsanyi), which in its original form did not take into account QCD massive splittings in the initial state. We constructed relevant kinematics and dipole splitting functions together with their integrals. We partially adapted the method in MC program and checked against the known massive inclusive NLO results.