Carbonate reservoirs typically exhibit very complex geological structures and are characterized by flow dynamics primarily occurring in fractures. The intricate network of fractures as well as their interconnectedness may lead to unexpected flow patterns and uneven sweep efficiency. Determining reservoir properties of both matrix and fracture channels is quintessential for accurately tracking the fluid front movement in the reservoir, optimizing sweep efficiency, and maximizing hydrocarbon production. In this study, we showcase the application of a feature-oriented ensemble-based history matching workflow to a complex fractured carbonate reservoir box model, focusing on the use of formation resistivity tomography data that are usually inferred from deep crosswell electromagnetic (EM) surveys. Compared with the production data that are commonly used in history matching, deep EM measurements provide additional information about the spatial distribution of subsurface reservoir properties in the interwell volumes by exploiting the strong resistivity contrast between water and hydrocarbons. A hybrid parameterization approach is used to represent the multiscale fracture distribution in which the spatial distribution of small-scale fractures is modelled by a truncated Gaussian simulation method. A large number (over one million) of uncertain model parameters including reservoir matrix and fracture properties as well as Archie's parameters are identified and updated by an iterative ensemble smoother. For an efficient integration of the high-dimensional and noisy EM tomography data, the boundary or contour information extracted from the EM resistivity field is instead assimilated through a distance parameterization approach. A modified bootstrap-based localization is proposed to regularize the model updates adaptively during the iteration to reduce sampling errors. Especially, to improve the computational efficiency in dealing with the large dimensions of both data and model parameters, the localization is implemented in a projected low-dimensional data subspace. Experimental results demonstrate the applicability and efficiency of the developed workflow for reservoir history matching in more realistic model settings. The comparative case study also illustrates the significance of jointly incorporating multiple sources of data for better quantification of model uncertainty, and the great potential of deep EM data for enhancing the characterization of complex fractured carbonate reservoirs.