Radiofrequency ablation (RFA) offers localized and minimally invasive thermal ablation of small-to-medium sized inoperable superficial tumors. In RFA, needles are inserted into the target with image guidance and current is passed through the needles, resulting in high temperatures and consequently target ablation. However, RFA has a high local recurrence rate caused by incomplete ablation. We therefore develop a novel two-stage mathematical framework for pre-operative inverse treatment planning where first, we identify needle positions and orientations using convex and integer programming techniques, referred to as needle orientation optimization (NOO), and then we determine the optimal thermal dose delivery for full target coverage by computing simultaneous thermal and electrostatic partial differential equations, referred to as thermal dose optimization (TDO). In NOO, different methodologies using geometric approximations for needle placement with and without trajectory considerations for single, multiple, and clustered RFA applicators are presented. Using outputs from NOO, in TDO, we perform thermal dose analysis, using several thermal damage models, for targets and organs-at-risk (OARs). Finally, we present scenario-based thermal damage analysis to understand the effect of translational and rotational needle deflection on target and OAR coverage. We test our framework on three clinical cases with four different margins, for a total of 12 targets. Our methodologies provide fast treatment plans that meet clinical guidelines, and our deflection analysis indicates that, depending on thermal damage model used, uncertain needle placement may significantly impede target coverage, and therefore, clinical study into causes of deflection are recommended.