How can I convert my matlab code into python?

close all;
clear all;
clc;

Datafiles = fileDatastore("temp_summary.12*.txt","ReadFcn",@readMonth,"UniformRead",true);
dataAll = readall(Datafiles);
dataAll.Year = year(dataAll.Day);
dataAll.Month = month(dataAll.Day);
dataAll.DD = day(dataAll.Day);
LY = (dataAll.Month(:)==2 & dataAll.DD(:)==29);
dataAll(LY,:) = [];

% Unstack variables
minT_tbl = unstack(dataAll,"MinT","Year","GroupingVariables", ["Month","DD"],"VariableNamingRule","preserve")
maxT_tbl = unstack(dataAll,"MaxT","Year","GroupingVariables", ["Month","DD"],"VariableNamingRule","preserve")
yrs =str2double(minT_tbl.Properties.VariableNames(3:end))';
% find min
[Tmin,idxMn] = min(minT_tbl{:,3:end},[],2,'omitnan');
Tmin_yr = yrs(idxMn);
% find max
[Tmax,idxMx] = max(maxT_tbl{:,3:end},[],2,'omitnan');
Tmax_yr = yrs(idxMx);
% find low high
[lowTMax,idxMx] = min(maxT_tbl{:,3:end},[],2,'omitnan');
LowTMax_yr = yrs(idxMx);
% find high low
[highlowTMn,idxMn] = max(minT_tbl{:,3:end},[],2,'omitnan');
HighLowT_yr = yrs(idxMn);

% find avg high
AvgTMx = round(mean(table2array(maxT_tbl(:,3:end)),2,'omitnan'));
% find avg low
AvgTMn = round(mean(table2array(minT_tbl(:,3:end)),2,'omitnan'));

% Results
tempTbl = [maxT_tbl(:,["Month","DD"]), table(Tmax,Tmax_yr,AvgTMx,lowTMax,LowTMax_yr,Tmin,Tmin_yr,AvgTMn,highlowTMn,HighLowT_yr)]
tempTbl2 = splitvars(tempTbl)
FID = fopen('Meda 12 Temperature Climatology.txt','w');
report_date = datetime('now','format','yyyy-MM-dd HH:MM');
fprintf(FID,'Meda 12 Temperature Climatology at %s \n', report_date);
fprintf(FID,"Month   DD   Temp Max (°F)  Tmax_yr   AvgTMax (°F)  lowTMax (°F)  LowTMax_yr   TempMin (°F)  TMin_yr   AvgTMin (°F)  HighlowTMin (°F)  HighlowT_yr \n");
fprintf(FID,'%3d %6d %7d %14d %11d %11d %15d %11d %13d %10d %13d %17d \n', tempTbl2{:,1:end}');
fclose(FID);
winopen('Meda 12 Temperature Climatology.txt')




function Tbl = readMonth(filename)
opts = detectImportOptions(filename);
opts.ConsecutiveDelimitersRule = 'join';
opts.MissingRule = 'omitvar';
opts = setvartype(opts,'double');
opts.VariableNames = ["Day","MaxT","MinT","AvgT"];
Tbl = readtable(filename,opts);
Tbl = standardizeMissing(Tbl,{999,'N/A'},"DataVariables",{'MaxT','MinT','AvgT'});
Tbl = standardizeMissing(Tbl,{-99,'N/A'},"DataVariables",{'MaxT','MinT','AvgT'});
[~,basename] = fileparts(filename);
% use the base file name, not the full file name:
d = str2double(extract(basename,digitsPattern));
if ~leapyear(d(3)) && d(2) == 2 % February of a non-leap-year
    Tbl(Tbl.Day == 29,:) = [];  % remove the 29th day data, if any
end
Tbl.Day = datetime(d(3),d(2),Tbl.Day);
end

function tf = leapyear(y)
if mod(y,4)       % year is not divisible by 4
    tf = false;   %   it is a common year
elseif mod(y,100) % year is not divisible by 100
    tf = true;    %   it is a leap year
elseif mod(y,400) % year is not divisible by 400
    tf = false;   %   it is a common year
else
    tf = true;    %   it is a leap year
end
end

[Walter Roberson]

... with difficulty.

[Jonathon Klepatzki]

This is what I have thus far....the issue that I am facing now is the last function using delimiter's rule and etc.

from glob2 import glob
import pandas as pd
from pandas import Series, DataFrame
from datetime import datetime
import matplotlib as plt
import sys
import time
import numpy as N
import math
# opening multiple file(s)
files = glob('*.txt','r') # this code works
files.year = files.year(files.day)
files.month = files.month(files.day)
files.dd = files.dd(files.dd)
LY = (files.month==2 & files.dd==29) #removing leap year data
files=N.array[LY,:] 
# unstack variables
minT_tbl=files.unstack("MinT","Year",)
maxT_tbl=files.unstack("MaxT","Year",)
yrs=float(minT_tbl.Properties.VariableNames[3,:])

# find min
[Tmin,idxMn]=N.amin(minT_tbl{:,3:end},[],2,'notnull')
Tmin_yr = yrs(idxMn)
# find max
[Tmax,idxMx]=N.amax(maxT_tbl{:,3:end},[],2,'notnull')
Tmax_yr = yrs(idxMx)
# find low high
[lowTmax,idxMx]=N.amin(maxT_tbl{:,3:end},[],2,'notnull')
LowTMax_yr = yrs(idxMx)
# find high low
[highlowTmn,idxMn]=N.amax(minT_tbl{:,3:end},[],2,'notnull')
HighLowT_yr = yrs(idxMn)

# find average high
AvgTMx=N.around(N.mean(pd.DataFrame(maxT_tbl{:,3:end})),2,'notnull')
# find average low
AvgTMn=N.around(N.mean(pd.DataFrame(minT_tbl{:,3:end})),2,'notnull')

#results
tempTbl=[maxT_tbl(:,["Month","dd"]),pd.DataFrame(Tmax,Tmax_yr,AvgTMx,lowTMax,LowTMax_yr,Tmin,Tmin_yr,AvgTMn,highlowTMn,HighLowT_yr)]
tempTbl2=N.split(tempTbl)
FID = open("Meda 12 Temperature Climatology.txt",'w')
report_date = time.ctime()
FID.write("Meda 12 Temperature Climatology at %s \n", report_date)
FID.write("Month   dd   Temp_Max (°F)  Tmax_yr   AvgTMax (°F)  lowTMax (°F)  LowTMax_yr   TempMin (°F)  TMin_yr   AvgTMin (°F)  HighlowTMin (°F)  HighlowT_yr \n")
FID.write("%3d %6d %7d %14d %11d %11d %15d %11d %13d %10d %13d %17d \n",tempTbl2{:,1:end}')
FID.close()
open("Meda 12 Temperature Climatology.txt")

def Tbl(filename):
opts = 

if ~leapyear(d(3)) & d(2) == 2:
Tbl(Tbl.Day == 29,:) = []
Tbl.Day = datetime(d(3),d(2),Tbl.Day)
def tf = leapyear(y):
if mod(y,4)
tf = False
elif mod(y,100)
tf = True
elif mod(y,400)
tf = False
else
tf = True

Answer 1

[surya venu]

Hi,

Regrettably, there isn't a straightforward method for converting MATLAB code directly into Python; such a process must be carried out manually. However, it's curious to know why you're considering this conversion when MATLAB provides the capability to integrate Python code. This feature allows you to utilize Python within MATLAB itself.

To know more about MATLAB Python Interface, check out the link below:

https://www.mathworks.com/help/matlab/call-python-libraries.html

Hope it helps.

[Jonathon Klepatzki]

Hi Surya,

The reason why I have to pursue Python is due to budgeting. The 2nd reason is because I am trying to learn this language for graduate school purposes.

Answer 2

[Carlos Dias Caisa]

% Flujo de carga trifasico y comparacion con el lineal
clc
clear all
% Tomar los datos
Feeder = LoadFeeder('FEEDER IEEE13_SINREGULADOR.xlsx');
% Dibujar
figure(1)
PlotFeeder(Feeder);
% Flujo de carga trifasico
Res = ThreePhase_LoadFlow(Feeder);
ShowResults(Res,Feeder);
function Feeder = LoadFeeder(Archivo);
% Load data from an Execel file
Feeder.Options.Name = Archivo;
Feeder.Topology = xlsread(Archivo,'Topology');
Feeder.Configurations = xlsread(Archivo,'Configurations');
Feeder.Loads = xlsread(Archivo,'Loads');
Feeder.Regulators = xlsread(Archivo,'Regulators');
General = xlsread(Archivo,'General');
Feeder.Slack = General(1);
if length(General(:,1)>=10)
    Feeder.Vpu_slack_phase = General(5:7).*exp(j*General(8:10)*pi/180);
else
    Feeder.Vpu_slack_phase = exp(j*[0,-120,120]*pi/180);    
end
Feeder.Vnom = General(2);
SistInterna = General(3);
Feeder.Switches = xlsread(Archivo,'Switches');
Feeder.Graphic = xlsread(Archivo,'Graphic');
Feeder.NumN = 0;
Feeder.NumL = length(Feeder.Topology(:,1));
Feeder.NumC = length(Feeder.Loads(:,1));
Feeder.NumZ = 0;
Feeder.NumR = 0;
Feeder.NumS = 0;
if length(Feeder.Switches)>0
    Feeder.NumS = length(Feeder.Switches(:,1));
end
if length(Feeder.Regulators)>0
    Feeder.NumR = length(Feeder.Regulators(:,1));
end

% Cambiar el nombre de los nodos
if Feeder.NumS>0
    Feeder.Nodes_ID = unique([Feeder.Topology(:,1);Feeder.Topology(:,2);Feeder.Switches(:,1);Feeder.Switches(:,2)]);
else
    Feeder.Nodes_ID = unique([Feeder.Topology(:,1);Feeder.Topology(:,2)]);    
end

Tpl = Feeder.Topology;
for k = 1:Feeder.NumL
    n = Tpl(k,1);
    nn = find(Feeder.Nodes_ID==n);
    Feeder.Topology(k,1) = nn;
    
    n = Tpl(k,2);
    nn = find(Feeder.Nodes_ID==n);
    Feeder.Topology(k,2) = nn;        
end

Swtc = Feeder.Switches;
for k = 1:Feeder.NumS
    n = Swtc(k,1);
    nn = find(Feeder.Nodes_ID==n);
    Feeder.Switches(k,1) = nn;
    
    n = Swtc(k,2);
    nn = find(Feeder.Nodes_ID==n);
    Feeder.Switches(k,2) = nn;    
end

% incluir los switches cerrados en las lineas
if Feeder.NumS>0
    sc = find(Feeder.Switches(:,3)==1);
    ls = length(sc);
    NM = zeros(ls,4);
    NM(:,1) = Feeder.Switches(sc,1);
    NM(:,2) = Feeder.Switches(sc,2);
    NM(:,3) = 0;  % distancia cero
    NM(:,4) = 1;  % configuracion 1 (no importa realmente porque la distancia es cero !
    Feeder.Topology = [Feeder.Topology; NM];
    Feeder.NumL = Feeder.NumL + ls;
end

NuevasCargas = zeros(Feeder.NumC,1);
for k = 1:Feeder.NumC
    n = Feeder.Loads(k,1);
    nn = find(Feeder.Nodes_ID==n);
    NuevasCargas(k) = nn;
end
Feeder.Loads(:,1) = NuevasCargas;
Feeder.Loads(:,4:9) = Feeder.Loads(:,4:9)*1000;  % pasar de kW a W
% Cambiar al sistema internacional si es el caso
if SistInterna==0
    Feeder.Topology(:,3) =  Feeder.Topology(:,3)/3280.8399;  % longitud pasa de ft a km
    Feeder.Configurations(:,3:20) = Feeder.Configurations(:,3:20)/1.609344;
end

% Impedancias de las lineas
Feeder.NumZ = length(Feeder.Configurations(:,1));
ZLIN = zeros(3,3,Feeder.NumZ);
BLIN = zeros(3,3,Feeder.NumZ);
for k = 1:Feeder.NumZ       
    ZLIN(1,1,k) = Feeder.Configurations(k,3)+j*Feeder.Configurations(k,9);
    ZLIN(1,2,k) = Feeder.Configurations(k,4)+j*Feeder.Configurations(k,10);
    ZLIN(2,1,k) = ZLIN(1,2,k);
    ZLIN(1,3,k) = Feeder.Configurations(k,5)+j*Feeder.Configurations(k,11);
    ZLIN(3,1,k) = ZLIN(1,3,k);
    ZLIN(2,2,k) = Feeder.Configurations(k,6)+j*Feeder.Configurations(k,12);
    ZLIN(2,3,k) = Feeder.Configurations(k,7)+j*Feeder.Configurations(k,13);
    ZLIN(3,2,k) = ZLIN(2,3,k);
    ZLIN(3,3,k) = Feeder.Configurations(k,8)+j*Feeder.Configurations(k,14);                
    
    BLIN(1,1,k) = Feeder.Configurations(k,15);
    BLIN(1,2,k) = Feeder.Configurations(k,16);
    BLIN(2,1,k) = BLIN(1,2,k);
    BLIN(1,3,k) = Feeder.Configurations(k,17);
    BLIN(3,1,k) = BLIN(1,3,k);
    BLIN(2,2,k) = Feeder.Configurations(k,18);
    BLIN(2,3,k) = Feeder.Configurations(k,19);
    BLIN(3,2,k) = BLIN(2,3,k);
    BLIN(3,3,k) = Feeder.Configurations(k,20);                        
end
Feeder.ZLIN = ZLIN;
Feeder.BLIN = 1/2*BLIN*1E-6;

% Importante:  Revisar por que es diferente:  No deberia
Feeder.NumN = length(Feeder.Nodes_ID);
if not(Feeder.NumN==Feeder.NumL+1)
    disp('Warning:  The system may be not radial because NumN != NumL+1');
end
for k = 1:Feeder.NumN
    n = Feeder.Graphic(k,1);
    nn = find(Feeder.Nodes_ID==n);
    NuevoGrafico(k,1) = nn;
end
Feeder.Graphic(:,1) = NuevoGrafico;
Feeder.Slack = find(Feeder.Nodes_ID==Feeder.Slack);

% cambiar los reguladores
if Feeder.NumR>0
    Feeder.Regulators(:,2:4) = Feeder.Regulators(:,2:4)*5/8/100 + 1;
end
Feeder.Options.DeltaLoadFlow = General(4);

function Ybust = Three_Phase_Ybus(Feeder);
    % matriz Ybus trifasica ordenada: primero todos los nodos de la fase A,
    % luego todos los nodos de la fase B y finalmente todos los nodos de C
    NumN = Feeder.NumN;
    Ybust = zeros(NumN*3);
    for k = 1:Feeder.NumL
        N1 = Feeder.Topology(k,1);
        N2 = Feeder.Topology(k,2);
        long = Feeder.Topology(k,3);
        c = Feeder.Topology(k,4);
        Zkm = Feeder.ZLIN(:,:,c)*long;    
        Ykm = inversa(Zkm);
        Bf = j*Feeder.BLIN(:,:,c)*long;
        kN1 = [N1,N1+NumN,N1+2*NumN];
        kN2 = [N2,N2+NumN,N2+2*NumN];
        Ybust(kN1,kN1)=Ybust(kN1,kN1)+Ykm;
        Ybust(kN1,kN2)=Ybust(kN1,kN2)-Ykm;
        Ybust(kN2,kN1)=Ybust(kN2,kN1)-Ykm;
        Ybust(kN2,kN2)=Ybust(kN2,kN2)+Ykm;
        Ybust(kN1,kN1)=Ybust(kN1,kN1)+Bf;
        Ybust(kN2,kN2)=Ybust(kN2,kN2)+Bf;
    end
        
function Y = inversa(Z)
    % calcula la inversa eliminando los ceros
    Y = zeros(3);
    W = sign(sum(abs(Z)));
    s = find(W==1);
    Yp = inv(Z(s,s));
    Y(s,s) = Yp;
    Vfslack = Feeder.Vpu_slack_phase*Vbase_fase; 
else
    Vfslack = Vbase_fase*exp([0,-j*120*pi/180,j*120*pi/180]);  % slack 1 < 0?
end

Vs(1,:) = Vfslack(1);
Vs(2,:) = Vfslack(2);
Vs(3,:) = Vfslack(3);

MFL = [1 -1 0; 0 1 -1; -1 0 1];% Cambia de fase a linea linea
%% barrido iterativo: con tensiones de fase
while (err>1E-10)
    Is = zeros(3,Feeder.NumN);   
    Vantes = Vs;
    Sperdidas = 0;
    % corrientes en las cargas: las corrientes van hacia abajo     
    for k = 1:Feeder.NumC
        N = Feeder.Loads(k,1);           
        estrella = Feeder.Loads(k,2);
        alpha = Feeder.Loads(k,3);
        Sestrella = zeros(3,1);
        Sdelta    = zeros(3,1);
        if estrella == 1
            Sestrella(1) = Feeder.Loads(k,4) + j*Feeder.Loads(k,5);
            Sestrella(2) = Feeder.Loads(k,6) + j*Feeder.Loads(k,7);
            Sestrella(3) = Feeder.Loads(k,8) + j*Feeder.Loads(k,9);
        else
            Sdelta(1) = Feeder.Loads(k,4) + j*Feeder.Loads(k,5);
            Sdelta(2) = Feeder.Loads(k,6) + j*Feeder.Loads(k,7);
            Sdelta(3) = Feeder.Loads(k,8) + j*Feeder.Loads(k,9);                 
        end
        Iest = conj(Sestrella./Vs(:,N)).*abs(Vs(:,N)/Vbase_fase).^alpha;
        Vlin = MFL*Vs(:,N);
        Idel = conj(Sdelta./Vlin).*abs(Vlin/Vbase_lin).^alpha;
        Is(:,N) = Is(:,N) + Iest + MFL'*Idel; % MFL debe estar transpuesta
    end
    Inodal = Is;       
    %% efecto capacitivo de las lineas
    for k = 1:Feeder.NumL
        N1 = Feeder.Topology(k,1);
        N2 = Feeder.Topology(k,2);
        long = Feeder.Topology(k,3);
        c = Feeder.Topology(k,4);
        Bf = Feeder.BLIN(:,:,c)*long;      
        Is(:,N1) = Is(:,N1) + j*Bf*Vs(:,N1);
        Is(:,N2) = Is(:,N2) + j*Bf*Vs(:,N2);
    end
    
    %% barrido de corrientes
    for k = 1:Feeder.NumL
        lin = Sort(k);
        N1 = Feeder.Topology(lin,1);
        N2 = Feeder.Topology(lin,2);
        r = LinReg(lin);        
        if  r== 0  % no es un regulador de tension
            Is(:,N1) = Is(:,N1) + Is(:,N2);
        else
            % el regulador esta en el envio          
            tp = Feeder.Regulators(r,2:4)';
            Is(:,N1) = Is(:,N1) + Is(:,N2)./tp;           
        end
    end
    %% barrido de tensiones
    for k = Feeder.NumL:-1:1
        lin = Sort(k);
        N1 = Feeder.Topology(lin,1);
        N2 = Feeder.Topology(lin,2);
        long = Feeder.Topology(lin,3);
        c = Feeder.Topology(lin,4);
        Zkm = Feeder.ZLIN(:,:,c)*long;
        Vkm = Zkm*Is(:,N2);
        r = LinReg(lin);
        if r > 0 % es un regulador            
            tp = Feeder.Regulators(r,2:4)';
            Vs(:,N2) = Vs(:,N1).*tp - Vkm;            
        else                
            Vs(:,N2) = Vs(:,N1) - Vkm;
        end       
        Sperdidas = Sperdidas+Is(:,N2)'*Vkm;
    end
    err = max(max(abs(Vs-Vantes)))/Vbase_fase;
    iter = iter + 1; 
    if iter > 100
        disp('Error: 100 iteraciones');
        break
    end
end
%% resultados
Vlin = conj(MFL*Vs)';

% Hacer cero las tensiones de las fases que no existen
for k = 1:Feeder.NumL           
    N2 = Feeder.Topology(k,2);
    c  = Feeder.Topology(k,4);    
    Zkm = abs(Feeder.ZLIN(:,:,c));
    H = sign(sum(Zkm))';    
    Vs(:,N2) = Vs(:,N2).*H;
    HF = 1-abs(MFL*H);   
    Vlin(N2,:) = Vlin(N2,:).*HF';
end

Res.Vpu_phase = conj(Vs')/Vbase_fase;
Res.Vpu_line  = Vlin/Vbase_lin;
Res.Inodes = -conj(Inodal');  % entrando al nodo
Res.Ilines =  conj(Is');
Res.iter = iter;
Res.err = err;
Res.Perd = real(Sperdidas/1000);
Res.Sort = Sort;
Res.iter = iter;

function Res = ThreePhase_LoadFlow(Feeder);
    % Unbalanced load flow in power distribution systems
    if Feeder.Options.DeltaLoadFlow==0
        Res = ThreePhase_LoadFlow_linetoneutral(Feeder); 
    else
        Feeder.Vslack_linetoline = [1 -1 0; 0 1 -1; -1 0 1]*Feeder.Vpu_slack_phase/sqrt(3);
        Res = ThreePhase_LoadFlow_delta(Feeder);    
    end
    
function ShowResults(Res,Feeder);
    
    fprintf('---------------%s---------------\n',Feeder.Options.Name);
    fprintf('Total Losses : %5.4f  kW\n',Res.Perd);
    %fprintf('N Iterations : %i  \n',Res.iter);
    if Feeder.Options.DeltaLoadFlow
        %% Voltajes de linea
        disp('------- LINE TO LINE VOLTAGES ----------------------------------------------------------');
        format short
        fprintf('NODE\tVAB(pu)\t VAB(deg)\t VBC(pu)\tVBC(deg)\tVCA(pu)\tVCA(deg)\n'); 
        for k = 1:Feeder.NumN
            V = abs(Res.Vpu_line(k,:));
            A = angle(Res.Vpu_line(k,:))*180/pi;
            n = Feeder.Nodes_ID(k);
            if (n<10)
                fprintf('N%i  \t',Feeder.Nodes_ID(k));
            else
                if (n<100)
                    fprintf('N%i \t',Feeder.Nodes_ID(k));    
                else
                    fprintf('N%i\t',Feeder.Nodes_ID(k));
                end
            end
            for m = 1:3     
                if (V(m)>0.1)
                    fprintf('%5.4f < ',V(m));
                    if (A(m)>0)
                        fprintf(' %5.4f \t ',A(m));
                    else
                        fprintf('%5.4f \t ',A(m));
                    end
                else
                    fprintf('.                  \t ') 
                end
            end
            fprintf('\n');            
        end
    else 
        % Voltajes de fase
        disp('--------PHASE VOLTAGES ---------------------------------------------------------');
        format short
        fprintf('NODE\tVAn(pu)\t VAn(deg)\t VBn(pu)\tVBn(deg)\tVCn(pu)\tVCn(deg)\n'); 
        for k = 1:Feeder.NumN
            V = abs(Res.Vpu_phase(k,:));
            A = angle(Res.Vpu_phase(k,:))*180/pi;
            n = Feeder.Nodes_ID(k);
            if (n<10)
                fprintf('N%i  \t',Feeder.Nodes_ID(k));
            else
                if (n<100)
                    fprintf('N%i \t',Feeder.Nodes_ID(k));    
                else
                    fprintf('N%i\t',Feeder.Nodes_ID(k));
                end
            end
            for m = 1:3     
                if (V(m)>0)
                    fprintf('%5.4f < ',V(m));
                    if (A(m)>0)
                        fprintf(' %5.4f \t ',A(m));
                    else
                        fprintf('%5.4f \t ',A(m));
                    end
                else
                    fprintf('.                  \t ') 
                end
            end
            fprintf('\n');            
        end
    end
    
function  PlotFeeder(Feeder);
    % Hacer el grafico
    NumN = length(Feeder.Graphic(:,1));
    NumL = length(Feeder.Topology(:,1));
    r = size(Feeder.Regulators);
    NumR = r(1);
    plot(Feeder.Graphic(:,2),-Feeder.Graphic(:,3),'.b');
    axis off
    hold on
    for k = 1:NumN
        text(Feeder.Graphic(k,2),-Feeder.Graphic(k,3),num2str(Feeder.Nodes_ID(Feeder.Graphic(k,1))));    
    end
    % Lineas
    for k = 1:NumL
        N1 = Feeder.Topology(k,1);
        N2 = Feeder.Topology(k,2);
        c  = Feeder.Topology(k,4);    
        Zkm = abs(Feeder.ZLIN(:,:,c));
        H = sum(sign(sum(Zkm))');        
        P1 = find(Feeder.Graphic(:,1)==N1);
        P2 = find(Feeder.Graphic(:,1)==N2);
        S1 =  Feeder.Graphic([P1,P2],2);
        S2 = -Feeder.Graphic([P1,P2],3);      
        if H==1
            plot(S1,S2,'m')    
        end
        if H==2
            plot(S1,S2,'r')    
        end
        if H==3
            plot(S1,S2,'b')    
        end
        
    end
    % Reguladores
    for k = 1:NumR
        R = Feeder.Regulators(k,1);    
        N1 = Feeder.Topology(R,1);
        N2 = Feeder.Topology(R,2);
        P1 = find(Feeder.Graphic(:,1)==N1);
        P2 = find(Feeder.Graphic(:,1)==N2);    
        SX = 0.5*(Feeder.Graphic(P1,2)+Feeder.Graphic(P2,2));
        SY = -0.5*(Feeder.Graphic(P1,3)+Feeder.Graphic(P2,3));
        plot(SX,SY,'ob','MarkerSize',12,'MarkerFaceColor',[0.7,0.7,0.9]);        
        text(SX-2,SY,'R');
    end
    % Switches
    if length(Feeder.Switches)>0
        for k = 1:length(Feeder.Switches(:,1))    
            N1 = Feeder.Switches(k,1);
            N2 = Feeder.Switches(k,2);
            P1 = find(Feeder.Graphic(:,1)==N1);
            P2 = find(Feeder.Graphic(:,1)==N2);    
            S1 = Feeder.Graphic([P1,P2],2);
            S2 = -Feeder.Graphic([P1,P2],3);
            SX = 0.5*(Feeder.Graphic(P1,2)+Feeder.Graphic(P2,2));
            SY = -0.5*(Feeder.Graphic(P1,3)+Feeder.Graphic(P2,3));            
            p = Feeder.Switches(k,3);
            if (p == 0) % switch abierto
                line(S1,S2,'Color',[.1 .8 .2]);            
                plot(SX,SY,'sg','MarkerFaceColor',[1,1,1]);        
            else
                line(S1,S2,'Color',[0 0 1]);
                plot(SX,SY,'sb');        
            end
        end    
    end
    title (Feeder.Options.Name)
    hold off

[Walter Roberson]

I do not understand this Answer? The task was to convert MATLAB code into Python, but the above posted code is in MATLAB, not in Python? The Answer also does not describe any procedure to convert MATLAB to Python ?