1010
1111import pysensors as ps
1212from matplotlib .patches import Circle
13- from pydmd import DMD
1413
1514
1615class GQR (QR ):
@@ -30,7 +29,7 @@ class GQR(QR):
3029
3130 @ authors: Niharika Karnik (@nkarnik2999), Mohammad Abdo (@Jimmy-INL), and Krithika Manohar (@kmanohar)
3231 """
33- def __init__ (self ,idx_constrained ,n_sensors ,n_const_sensors ,all_sensors ,constraint_option ,nx ,ny ,r ):
32+ def __init__ (self ): # ,idx_constrained,n_sensors,n_const_sensors,all_sensors,constraint_option,nx,ny,r
3433 """
3534 Attributes
3635 ----------
@@ -49,21 +48,18 @@ def __init__(self,idx_constrained,n_sensors,n_const_sensors,all_sensors,constrai
4948 exact_n_const_sensors : The number of sensors in the constrained region should be exactly equal to n_const_sensors.
5049 """
5150 self .pivots_ = None
52- self .optimality = None
53-
54- self .constrainedIndices = idx_constrained
55- self .nSensors = n_sensors
56- self .nConstrainedSensors = n_const_sensors
57- self .all_sensorloc = all_sensors
58- self .constraint_option = constraint_option
59- self ._nx = nx
60- self ._ny = ny
61- self ._r = r
62-
63- def fit (
64- self ,
65- basis_matrix
66- ):
51+ # self.optimality = None
52+
53+ # self.constrainedIndices = idx_constrained
54+ # self.n_sensors = n_sensors
55+ # self.nConstrainedSensors = n_const_sensors
56+ # self.all_sensorloc = all_sensors
57+ # self.constraint_option = constraint_option
58+ # self._nx = nx
59+ # self._ny = ny
60+ # self._r = r
61+
62+ def fit (self ,basis_matrix = None ,** optimizer_kws ):
6763 """
6864 Parameters
6965 ----------
@@ -77,16 +73,48 @@ def fit(
7773 -------
7874 self: a fitted :class:`pysensors.optimizers.QR` instance
7975 """
76+ if 'idx_constrained' in optimizer_kws .keys ():
77+ self .constrainedIndices = optimizer_kws ['idx_constrained' ]
78+ else :
79+ self .constrainedIndices = []
80+ if 'n_sensors' in optimizer_kws .keys ():
81+ self .n_sensors = optimizer_kws ['n_sensors' ]
82+ else :
83+ self .n_sensors = np .shape (basis_matrix )[0 ]
84+ if 'n_const_sensors' in optimizer_kws .keys ():
85+ self .nConstrainedSensors = optimizer_kws ['n_const_sensors' ]
86+ else :
87+ self .nConstrainedSensors = None
88+ if 'all_sensors' in optimizer_kws .keys ():
89+ self .all_sensors = optimizer_kws ['all_sensors' ]
90+ else :
91+ self .all_sensors = None
92+ if 'constraint_option' in optimizer_kws .keys ():
93+ self .constraint_option = optimizer_kws ['constraint_option' ]
94+ else :
95+ self .constraint_option = None
96+ if 'nx' in optimizer_kws .keys ():
97+ self ._nx = optimizer_kws ['nx' ]
98+ else :
99+ self ._nx = None
100+ if 'ny' in optimizer_kws .keys ():
101+ self ._ny = optimizer_kws ['ny' ]
102+ else :
103+ self ._ny = None
104+ if 'r' in optimizer_kws .keys ():
105+ self ._r = optimizer_kws ['r' ]
106+ else :
107+ self ._r = None
80108
81109 n_features , n_samples = basis_matrix .shape # We transpose basis_matrix below
82- max_const_sensors = len (self .constrainedIndices ) #Maximum number of sensors allowed in the constrained region
110+ max_const_sensors = len (self .constrainedIndices ) # Maximum number of sensors allowed in the constrained region
83111
84112 ## Assertions and checks:
85- # if self.nSensors > n_features - max_const_sensors + self.nConstrainedSensors:
113+ # if self.n_sensors > n_features - max_const_sensors + self.nConstrainedSensors:
86114 # raise IOError ("n_sensors cannot be larger than n_features - all possible locations in the constrained area + allowed constrained sensors")
87- # if self.nSensors > n_samples + self.nConstrainedSensors: ## Handling zero constraint?
115+ # if self.n_sensors > n_samples + self.nConstrainedSensors: ## Handling zero constraint?
88116 # raise IOError ("Currently n_sensors should be less than min(number of samples, number of modes) + number of constrained sensors,\
89- # got: n_sensors = {}, n_samples + const_sensors = {} + {} = {}".format(self.nSensors ,n_samples,self.nConstrainedSensors,n_samples+self.nConstrainedSensors))
117+ # got: n_sensors = {}, n_samples + const_sensors = {} + {} = {}".format(self.n_sensors ,n_samples,self.nConstrainedSensors,n_samples+self.nConstrainedSensors))
90118
91119 # Initialize helper variables
92120 R = basis_matrix .conj ().T .copy ()
@@ -96,33 +124,38 @@ def fit(
96124
97125 for j in range (k ):
98126 r = R [j :, j :]
127+
99128 # Norm of each column
100129 dlens = np .sqrt (np .sum (np .abs (r ) ** 2 , axis = 0 ))
130+
101131 if self .constraint_option == "max_n_const_sensors" :
102- dlens_updated = ps .utils ._norm_calc .norm_calc_max_n_const_sensors (self .constrainedIndices ,dlens ,p ,j , self .nConstrainedSensors ,self .all_sensorloc ,self .nSensors )
132+ dlens_updated = ps .utils ._norm_calc .norm_calc_max_n_const_sensors (self .constrainedIndices ,dlens ,p ,j , self .nConstrainedSensors ,self .all_sensorloc ,self .n_sensors )
103133 i_piv = np .argmax (dlens_updated )
104134 dlen = dlens_updated [i_piv ]
105- elif self .constraint_option == "exact_n_const_sensors" :
135+ elif self .constraint_option == "exact_n_const_sensors" :
106136 dlens_updated = ps .utils ._norm_calc .norm_calc_exact_n_const_sensors (self .constrainedIndices ,dlens ,p ,j ,self .nConstrainedSensors )
107137 i_piv = np .argmax (dlens_updated )
108138 dlen = dlens_updated [i_piv ]
109139 elif self .constraint_option == "predetermined_end" :
110- dlens_updated = ps .utils ._norm_calc .predetermined_norm_calc (self .constrainedIndices , dlens , p , j , self .nConstrainedSensors , self .nSensors )
140+ dlens_updated = ps .utils ._norm_calc .predetermined_norm_calc (self .constrainedIndices , dlens , p , j , self .nConstrainedSensors , self .n_sensors )
111141 i_piv = np .argmax (dlens_updated )
112142 dlen = dlens_updated [i_piv ]
113143 elif self .constraint_option == "radii_constraints" :
114-
115- if j == 0 :
144+
145+ if j == 0 :
116146 i_piv = np .argmax (dlens )
117147 dlen = dlens [i_piv ]
118148 dlens_old = dlens
119149 else :
120150
121- dlens_updated = ps .utils ._norm_calc .f_radii_constraint (j ,dlens ,dlens_old ,p ,self ._nx ,self ._ny ,self ._r ,self .all_sensorloc , self .nSensors ) #( self.radius,self._nx,self._ny,self.all_sensorloc,dlens,p,j)
151+ dlens_updated = ps .utils ._norm_calc .f_radii_constraint (j ,dlens ,dlens_old ,p ,self ._nx ,self ._ny ,self ._r ,self .all_sensorloc , self .n_sensors ) #( self.radius,self._nx,self._ny,self.all_sensorloc,dlens,p,j)
122152 i_piv = np .argmax (dlens_updated )
123153 dlen = dlens_updated [i_piv ]
124154 dlens_old = dlens_updated
125-
155+ else :
156+ i_piv = np .argmax (dlens )
157+ dlen = dlens [i_piv ]
158+
126159 # Choose pivot
127160 # i_piv = np.argmax(dlens_updated)
128161
@@ -148,7 +181,6 @@ def fit(
148181 R [j + 1 :, j ] = 0
149182
150183 self .pivots_ = p
151-
152184 return self
153185
154186if __name__ == '__main__' :
@@ -191,10 +223,10 @@ def plot_gallery(title, images, n_col=n_col, n_row=n_row, cmap=plt.cm.gray):
191223
192224 #Find all sensor locations using built in QR optimizer
193225 #max_const_sensors = 230
194- n_const_sensors = 3
195- n_sensors = 4
196- n_modes = 40
197- r = 2
226+ n_const_sensors = 0
227+ n_sensors = 10
228+ n_modes = 10
229+ r = 5
198230 # dmd = DMD(svd_rank=0,exact=True,opt=False)
199231 # dmd.fit(X.T)
200232 # U = dmd.modes.real
@@ -209,9 +241,9 @@ def plot_gallery(title, images, n_col=n_col, n_row=n_row, cmap=plt.cm.gray):
209241 # top_sensors_dmd_unconstrained = model_dmd_unconstrained.get_selected_sensors()
210242 # optimality_dmd = ps.utils._validation.determinant(top_sensors_dmd_unconstrained, n_features, basis_matrix_dmd)
211243 # print(optimality0)
212- # basis = ps.basis.SVD(n_basis_modes=n_modes)
244+ basis = ps .basis .SVD (n_basis_modes = n_modes )
213245 optimizer = ps .optimizers .QR ()
214- model = ps .SSPOR (optimizer = optimizer , n_sensors = n_sensors )
246+ model = ps .SSPOR (optimizer = optimizer , n_sensors = n_sensors , basis = basis )
215247 model .fit (X_small )
216248 top_sensors0 = model .get_selected_sensors ()
217249 all_sensors = model .get_all_sensors ()
@@ -227,32 +259,32 @@ def plot_gallery(title, images, n_col=n_col, n_row=n_row, cmap=plt.cm.gray):
227259 # didx = np.isin(all_sensors,sensors_constrained,invert=False)
228260 # const_index = np.nonzero(didx)
229261 # j =
230-
231262
232263
233- ##Plotting the constrained region
234- # ax = plt.subplot()
235- # #Plot constrained space
236- # img = np.zeros(n_features)
237- # img[sensors_constrained] = 1
238- # im = plt.imshow(img.reshape(image_shape),cmap=plt.cm.binary)
239- # # create an axes on the right side of ax. The width of cax will be 5%
240- # # of ax and the padding between cax and ax will be fixed at 0.05 inch.
241- # divider = make_axes_locatable(ax)
242- # cax = divider.append_axes("right", size="5%", pad=0.05)
243- # plt.colorbar(im, cax=cax)
244- # plt.title('Constrained region');
264+
265+ #Plotting the constrained region
266+ ax = plt .subplot ()
267+ #Plot constrained space
268+ img = np .zeros (n_features )
269+ img [sensors_constrained ] = 1
270+ im = plt .imshow (img .reshape (image_shape ),cmap = plt .cm .binary )
271+ # create an axes on the right side of ax. The width of cax will be 5%
272+ # of ax and the padding between cax and ax will be fixed at 0.05 inch.
273+ divider = make_axes_locatable (ax )
274+ cax = divider .append_axes ("right" , size = "5%" , pad = 0.05 )
275+ plt .colorbar (im , cax = cax )
276+ plt .title ('Constrained region' )
245277
246278 ## Fit the dataset with the optimizer GQR
247279 optimizer1 = GQR (sensors_constrained ,n_sensors ,n_const_sensors ,all_sensors , constraint_option = "radii_constraints" ,nx = nx , ny = ny , r = r )
248- model1 = ps .SSPOR ( optimizer = optimizer1 , n_sensors = n_sensors )
280+ model1 = ps .SSPOR ( optimizer = optimizer1 , n_sensors = n_sensors , basis = basis )
249281 model1 .fit (X_small )
250282 all_sensors1 = model1 .get_all_sensors ()
251283 basis_matrix = model1 .basis_matrix_
252284 top_sensors = model1 .get_selected_sensors ()
253285 print (top_sensors )
254- # optimality = ps.utils._validation.determinant(top_sensors, n_features, basis_matrix)
255- # print(optimality)
286+ optimality = ps .utils ._validation .determinant (top_sensors , n_features , basis_matrix )
287+ print (optimality )
256288 # optimizer_dmd_constrained = ps.optimizers.GQR(sensors_constrained,n_sensors,n_const_sensors,all_sensors_dmd_unconstrained,constraint_option = "exact_n_const_sensors",nx = nx, ny = ny, r = r)
257289 # model_dmd_constrained = ps.SSPOR(n_sensors=n_sensors, basis=ps.basis.Custom(n_basis_modes=n_modes, U=U), optimizer = optimizer_dmd_constrained)
258290 # model_dmd_constrained.fit(X)
@@ -281,14 +313,21 @@ def plot_gallery(title, images, n_col=n_col, n_row=n_row, cmap=plt.cm.gray):
281313 img = np .zeros (n_features )
282314 img [top_sensors ] = 16
283315 fig ,ax = plt .subplots (1 )
284- ax . set_aspect ( 'equal' )
316+
285317 ax .imshow (img .reshape (image_shape ),cmap = plt .cm .binary )
286318 print (top_sensors )
287319 top_sensors_grid = np .unravel_index (top_sensors , (nx ,ny ))
288320 # figure, axes = plt.subplots()
289321 for i in range (len (top_sensors_grid [0 ])):
290322 circ = Circle ( (top_sensors_grid [1 ][i ], top_sensors_grid [0 ][i ]), r ,color = 'r' ,fill = False )
291323 ax .add_patch (circ )
324+ # ax.plot([xmin,xmin],[ymin,ymax],'-r')
325+ # ax.plot([xmax,xmax],[ymin,ymax],'-r')
326+ # ax.plot([xmin,xmax],[ymin,ymin],'-r')
327+ # ax.plot([xmin,xmax],[ymax,ymax],'-r')
328+ ax .set_aspect ('equal' )
329+ # ax.set_xlim([0,64])
330+ # ax.set_ylim([0,64])
292331 plt .show ()
293-
294-
332+
333+
0 commit comments