So, an important thing to know here is that a normal is only a vector, so it, in a way, represents a single axis of the rotation. The script below pretty much guesses an x and y heading and uses the z direction as the up axis for the cursor. Here, I'm finding the location and normal of the selected face and assigning it to the location and up axis of the cursor. I've simplified this script a lot. It doesn't use the bmesh data, which means you may need to make your face selection in edit mode, and switch back to object mode and then back into edit mode before you run the script.
import bpy
import mathutils
C = bpy.context
cursor = C.scene.cursor
obj = C.active_object
p = obj.data.polygons
#get the normal and center of the selected face
normal = obj.data.polygons[p.active].normal
center = obj.data.polygons[p.active].center
#generate a quaternion rotation from the normal vector
up = normal.to_track_quat('X', 'Y')
#set the cursor to the center location and rotate it to match the up vector
cursor.location = center
cursor.rotation_quaternion = up
And then back in object mode again carry on with getting the cursor location and rotation and creating the primative like in your example above:
#get the object world matrix
mat = obj.matrix_world
#generate the location and rotation
rot = cursor.rotation_euler.copy()
loc = mat @ cursor.location.copy()
bpy.ops.object.mode_set(mode='OBJECT')
bpy.ops.mesh.primitive_plane_add(location = loc, rotation = rot)
You have to get the object world matrix and use it to offset the center location because any vertex, edge, or face location will be based on a distance from the object center, not the world center. Of course, you don't have to use the cursor as a step in between and can just apply the rotation and location of the face directly using:
rot = up
loc = mat @ center