Portfolio optimization with a conditional value in Python

I have the covariance matrix, the return vector and some scores (ESG scores). The objective is to derive the efficient frontier, the minimum variance portfolio and the tangency portfolio conditioned on the average ESG score of the portfolio being above a certain treshold. I am able to produce the usual efficient frontier but I am not sure how to implement the condition of the ESG score. Most of the code is provided by The Quant MBA blog here.

Below is my code in Python. Note that the variable avg_esg is what the average score of the tangency portfolio is found to be, but here I have not set any condition beforehand so this is merely what it ends up with optimizing without the condition. In my original code this is wrapped in a function that calls a dataframe, but for reproduceability I have just provided the matrices of a sample I am testing on.

# imports
from __future__ import division
from matplotlib import pyplot as plt
from numpy.linalg import inv,pinv
import numpy as np

""" USING THE MATRIX NOTATION DEFINED IN
APPENDIX OF "A CRITIQUE OF THE ASSET PRICING THEORY'S TESTS" ROLL (1977) 
"""

# Source code can be found at: 
# https://thequantmba.wordpress.com/2016/04/16/modern-portfolio-theory-python/

# USER INPUT
# V covariance matrix, R is the return vector
V = np.matrix([[84.76695659,  20.8854772,   20.62182415,  74.73652696,  14.35995947], 
              [20.8854772,   35.22429277,  12.95439707,  32.22912903,  12.96449085],
              [20.62182415,  12.95439707,  44.02079739,  38.73627316,   9.46608475],
              [74.73652696,  32.22912903,  38.73627316, 178.86640813,  33.40281336],
              [14.35995947,  12.96449085,   9.46608475,  33.40281336,  32.38514103]])
R = np.matrix([[-0.32264539, -0.08469428, 1.27628749, -0.23207085, 0.21012106]]).T
rf =  0.1 # risk-free rate 

# Mean ESG score of each company
esgarr = np.matrix([[8.24336898, 4.6373262,  8.30657754, 4.65406417, 3.43620321]]).T

# define the vector of ones and standard deviations
I = np.ones((len(R),1)) # vector of ones 
SD = np.sqrt(np.diag(V)) # standard deviation of each company

# Variable for Efficient frontier (SEE (A.9) ROLL 1977 page 160)
# The efficient set constants contained in the information matrix
# The paper describes the notation 
# Link to free access paper: https://www.scinapse.io/papers/2011043606 
C = I.T*pinv(V)*I 
B = R.T*pinv(V)*I
A = R.T*pinv(V)*R 
# D is the determinant used for evaluating if it is the global min/max/saddle
D = A*C-B**2

########################################################
#EFFICIENT FRONTIER
########################################################
mu = np.arange(-max(R),max(R)*5,max(R)/100)
minvar = (A-2*B*mu+(C*mu**2))/D;
minstd = np.sqrt(minvar)[0];
minstd = np.squeeze(np.asarray(minstd))

########################################################
#MVP
########################################################
# Mean and Variance of Global Minimum Variance Portfolio
mu_g = B/C
var_g = 1/C
std_g = np.sqrt(var_g)

# Minimum Variance Portfolio Weights
w_g = (pinv(V)*I)/C

########################################################
#TANGENT PORTFOLIO
########################################################

# Expected Return of Tangency Portfolio
mu_tan = (B*rf-A)/(C*rf-B);

# Variance and Standard Deviation of Tangency Portfolio
vartan = (A-2*rf*B + (rf**2*C))/((B-C*rf)**2);
stdtan = np.sqrt(vartan);

# Weights for Tangency Portfolio
w_tan = (pinv(V)*(R - rf*I))/(B-C*rf)
# ESG weights and average esg score
w_esg = np.dot(w_tan.T, esgarr)
w_tan_i = np.dot(w_tan.T, I)
avg_esg = np.divide(w_esg, w_tan_i)
print('Weights of each security in the tangency portfolio \n', w_tan.T
      ,'\nSum of weights', sum(w_tan))
print('Average ESG score for each security \n', esgarr.T)

# Tangency Line
m_tan = mu[mu >= rf];
minvar_rf = (m_tan-rf)**2/(A-2*rf*B+C*rf**2);
minstd_rf = np.sqrt(minvar_rf);
# np.squeeze removes single-dimensional entries from shape of array
minstd_rf = np.squeeze(np.asarray(minstd_rf)) 
sharpe = mu_tan / stdtan

# Plot with tangency portfolio
plt.plot(minstd,mu,SD,R,'*',minstd_rf,m_tan,'r',std_g,mu_g,'bo',stdtan,mu_tan,'go')
plt.title('Efficient Frontier with Individual Securities',fontsize=18)
plt.ylabel('Expected Return (%)', fontsize=12)
plt.xlabel('Standard Deviation (%)', fontsize=12)
plt.text(1,rf,'rf',fontsize=12);
plt.text(1+std_g,mu_g,'Min-var portfolio',fontsize=12);
plt.text(1+stdtan,mu_tan,'Tangency Portfolio',fontsize=12);
#plt.savefig("eff_frontier.png")
plt.show()
print('Sharpe ratio of the tangency portfolio is', sharpe)
print('Average ESG score of portfolio is', avg_esg[0])
#print(m_tan)

from Portfolio optimization with a conditional value in Python