% VisNetDriver3.m in VisNetMat\VisNetMat3\ 9 December 2025
% The functions called are CreateNet.m, train.m, test.m, SetIMLIST.m, gui.m, and info.m (which calls infors_visnet.m and inform_visnet.m)
% The documentation is in VisNetMat3Doc25.pdf
% After a run has finished, type 'gui' in the command window to view the results
% In the gui, do not alter 'Loc' from its value of 1, for translation to different locations is not implemented and needs to be learned.
% For literature see Rolls,E.T. (2023) Brain Computations and Connectivity. Oxford University Press. Chapter 2 and Section D.6 and
% Zhang C, Rolls ET, Feng J. (2026) Invariant visual object and face learning in the ventral cortical visual pathway: a biologically plausible
% model. PLoS Computational Biology and
% Rolls,E.T. (2026) Brain Computations and Principles; and AI. Oxford University Press. Chapter 2 and Section D.6 and
% papers at https://www.oxcns.org
% copyright E T Rolls, Oxford Centre for Computational Neuroscience, https://www.oxcns.org
% No warranty express or implied is provided: this is research software.
% If this software is used or adapted for new published research, 
% it is requested that the original publications listed above, which provide details, are cited.
%
% % This program expects to find a directory VNanalysis at the same level
% as this program, so that ..\VNanalysis works.
% If this program has not been run before for a given NET_SIZE, then CreateNetFlag should be set to 1, 
% so that an Untrained net is created.
% Windows 11 requires \ in its paths, and Linux requires /, so edit as necessary below and in SetIMLIST.m.
% This program also expects to find base_dir defined correctly in SetIMLIST.m to find the imagedata directory.

clear; % this should leave compiled code ready for reuse
clear global; % make sure
close all hidden; format compact;
t1 = clock;

global EFR EFRDATA NLAYERS NET_SIZE SQNET_SIZE N_FILES N_GROUPS N_VIEWS N_LOCS BESTCELLS NBESTCELLS ...
       NBESTCELLSINFORM NTRANSFORMS PARAMETER_LIST SYNMAT SYNMAT_CHANGE ...
       SYNMAT_CHANGE_GROUP LOGNUM FIXEDTRACE SPARSENESS BETA ETA LRNRATE NSYN SYNMATL1 SYNMATL2 SYNMATL3 ...
       SYNMATL4 VNLOGINDEX IMLIST IMLIST_TEST N_VIEWS_TEST FILTIM NSYNL1 Images base_dir image_dir ...
       FiltImSize SqFiltImSize Display ...
       filter_layer_1 filter_layer_2  filter_layer_3 filter_layer_4 SIGMOID LATERALINHIB TraceDisp DisplaySYNMATprevious ...
       NORM MAX_WEIGHTS
   
NLAYERS = 4;        % number of layers of VisNet
NET_SIZE = 32;      % size of network: 16 or 32 or 64 or 128 or 256
SQNET_SIZE = NET_SIZE * NET_SIZE;
N_GROUPS = 9;       % number of objects. 3 for ImageSet=1 (defined in SetIMLIST.m); 9 or up to 50 for ImageSet=2
N_VIEWS = 9;        % number of transforms of each object for learning. Overwritten in SetIMLIST
N_VIEWS_TEST = 9;   % number of transforms for testing: 9 for objects, 3 for bars. 3 for ImageSet=1; 9 for ImageSet=2. Overwritten in SetIMLIST
Interpolation = 0;  % Set this below to 1 to test with interpolated images
SetIMLIST(IMLIST, IMLIST_TEST, base_dir, Interpolation); % Set up the file names for the training and test objects in SetIMLIST.m
N_LOCS = 1; % This version uses only this value
locs = 1;   % and this is needed
EFR = zeros(SQNET_SIZE, NLAYERS);       % firing rates for layers 1 - 4
EFRDATA = zeros(N_GROUPS, N_VIEWS_TEST, NLAYERS, SQNET_SIZE); % the data saved in SweepFilename and used for testing the network
NBESTCELLS = 5;                         % the number of best cells for each stimulus to use for the information analysis
BESTCELLS = zeros(N_GROUPS * NBESTCELLS, 3); % rows = 5 cells per stim; col1=cell  col 2=stim  col3=info
SPARSENESS = [0.01 0.01 0.01 0.01] * 1.0;% the sparseness for each layer. 0.01 works for NET_SIZE 32 to 256
if NET_SIZE == 16
    SPARSENESS = [0.015 0.015 0.015 0.015];     % the sparseness for each layer. Works for NET_SIZE 16 so that we train sufficient cells
end
ETA = [0.0 0.8 0.8 0.8];                % the trace value for each layer
LRNRATE = [0.1 0.1 0.1 0.1] * 0.5;      % for the synaptic modification rule (default 0.1). When preRates is scaled for layer 1 to 0-1
% NSYN = [272 200 200 200];             % synapses per neuron
% NSYNL1 = [201 50 13 8];               % the number of L1 synapses for the different frequency filters, highest first
NSYN = [340 200 200 200] * 1;           % synapses per neuron for VisNetMat2
NSYNL1 = [128 32 8 2] * 2;              % the number of L1 synapses for the different frequency filters, highest first VisNetMat2
FiltImSize = 512;                       % the size of the filtered image files is 512x512
SqFiltImSize = FiltImSize * FiltImSize;
SIGMOID = 1;                            % set this to 1 for sigmoid activation function, 0 for linear threshold 
BETA = [10 10 10 10] * 1;               % beta, the slope of the sigmoid activation function for each layer. [1 1 1 1]
LATERALINHIB = 1;                       % 1 for convolution of the firing rates in a layer with a lateral inhibition filter. 0 for none.

% The following 2 are new for VisNet2
NORM = 3;                               % 1 for weight vector normalisation; 2 for Oja rule; 3 for Hebb rule with Long-Term Depression; 4 for MAX_WEIGHTS only
MAX_WEIGHTS = [0.1 0.1 0.1 0.1] * 0;    % if zero, do not set MAX_WEIGHTS on each neuron. else set them to these values for each layer
%MAX_WEIGHTS = [0.0 0.0 0.0 0.0]; 
 

% New nets with new random connections will be generated when CreateNetFlag=1
CreateNetFlag = 1;                      % set to 1 create a new net, or 0 to load a previously trained net. For NET_SIZE=16, leave this as 1, because it is fast.
TrainEpochs = [7 7 7 7] * 1;            % set to the number of epochs to train each of the 4 layers. At least [7 7 7 7] recommended for quantitative studies, to allow for convergence. No training if all 0.
disp(['Matlab will train with epochs for layers 1-4 = ' int2str(TrainEpochs)]);
% ImageSet = 2; % 1 is for bars; 2 is for ALOI objects. You can generate more. Set this in SetIMLIST.m
%rng('shuffle'); % set random num gen by the time of day. This is generally recommended, with the results of several runs being checked. 
rng(2); % reset the random number generator used in Matlab, only if you wish to create the same net
        % for example so that the same cells learn so that their trace can be viewed.
        
% Most of the parameters below this need not be changed, except perhaps for receptiveFieldWidth

Nfreq = 4;                              % the input images have been Gabor filtered to emulate V1
Nang = 4;                               % with this number of frequencies etc. See section 4 of VisNet_C_Matlab20.docx
Nsign = 2;
SYNMATL1 = spalloc(SqFiltImSize * Nfreq * Nang * Nsign, SQNET_SIZE, NSYN(1) * SQNET_SIZE); % A Matlab sparse matrix for the synaptic weights
SYNMATL2 = spalloc(SQNET_SIZE, SQNET_SIZE, NSYN(2) * SQNET_SIZE);
SYNMATL3 = spalloc(SQNET_SIZE, SQNET_SIZE, NSYN(3) * SQNET_SIZE);
SYNMATL4 = spalloc(SQNET_SIZE, SQNET_SIZE, NSYN(4) * SQNET_SIZE);

receptiveFieldWidth = [15 7 7 7];       % the width of the Gaussian region from which a neuron receives connections. It is scaled for FiltImSize. Was [15 7 7 7]
%receptiveFieldWidth = [15 12 12 12];       % the width of the Gaussian region from which a neuron receives connections. It is scaled for FiltImSize. Was [15 7 7 7]
receptiveFieldWidth = receptiveFieldWidth * (NET_SIZE / 32); % This scales for different net sizes so that the previous line need not be changed
receptiveFieldWidth(1) = 15;            % this does not change for different net sizes. It was 25, and that is still worth trying.

layer0 = [FiltImSize FiltImSize];
layer1 = [NET_SIZE NET_SIZE];
fieldWidth0 = floor(layer0 / 2);
fieldWidth = floor(layer1 / 2);

FILTIM = zeros(1, FiltImSize * FiltImSize* Nfreq * Nang * Nsign);  % 4 freqs, 4 angles, pos and neg

VNLOGINDEX = 1;                         % This can be manually updated to specify the names of files saved to create a record of experiments 
Display = 1;                            % 0=no display; 1=some display; 2=all display

TrainMode = 2;     % 1=translation invariance   2=translation invariance and view invariance. Only 2 is implemented.

% The lateral inhibition filter used if LATERALINHIB = 1. A sum 0 filter will produce balanced neg and pos value
load dog1_0.2_4.mat dog % difference of gaussian filter scale 1 centre sd=0.2 surround sd=4 size 31x31 for the lateralinhibition. See dogeg.m
filter_layer_1 = dog; % the new filters for use with Mode 2 and conv2 made by dogeg.m
filter_layer_2 = dog; % 
filter_layer_3 = dog;
filter_layer_4 = dog;

if NET_SIZE == 32
    VNNetLoadName ='..\VNanalysis\netUntrained32.mat'; % Beware: Windows uses \; Linux uses /
    VNNetTrainedLoadName = '..\VNanalysis\netTrained32.mat';
elseif NET_SIZE == 64
    VNNetLoadName ='..\VNanalysis\netUntrained64.mat';
    VNNetTrainedLoadName = '..\VNanalysis\netTrained64.mat';
elseif NET_SIZE == 128
    VNNetLoadName ='..\VNanalysis\netUntrained128.mat';
    VNNetTrainedLoadName = '..\VNanalysis\netTrained128.mat';
elseif NET_SIZE == 256
    VNNetLoadName ='..\VNanalysis\netUntrained256.mat';
    VNNetTrainedLoadName = '..\VNanalysis\netTrained256.mat';
else
    VNNetLoadName ='..\VNanalysis\netUntrained.mat';
    VNNetTrainedLoadName = '..\VNanalysis\netTrained.mat';
end
VNNetTrainedLoadName = 0; % comment this out if you want to load a trained net for testing
SYNMAT_CHANGE_GROUP = ones(NLAYERS, (max(TrainEpochs)) * N_GROUPS); % one every time that a new group is trained

NTRANSFORMS = N_VIEWS;  

if CreateNetFlag == 1       % Create new net. Otherwise load a previously created untrained net to save time
    CreateNet(receptiveFieldWidth, fieldWidth0, fieldWidth, Nfreq, Nang, Nsign); % CreateNet2 scales up the input to fill the size of VisNet
else
    load(char(VNNetLoadName));
    disp(['Net to be loaded is ', VNNetLoadName, '  NSYN =', num2str(NSYN)]);
end

disp(['VisNet log number= ' int2str(VNLOGINDEX)]);
LOGNUM = VNLOGINDEX;
VNNetSaveName = ['..\VNanalysis\net', int2str(VNLOGINDEX)];
SweepFilename = ['..\VNanalysis\Sweep', int2str(VNLOGINDEX), '.mat'];
% tmp = strcat('cp VisNetDriver7.m ../VNanalysis/VisNetDriver7_', int2str(VNLOGINDEX), '.m'); % linux
%tmp = strcat('copy VisNetDriver7b.m ..\VNanalysis\VisNetDriver7b_', int2str(VNLOGINDEX), '.m') % Windows 10
%system(tmp);

if sum(TrainEpochs) > 0    % if you do not wish to train, only test, then set all entries in TrainEpochs to 0
    
    time1 = clock;
    
    for layer = 1 : NLAYERS
        if TrainEpochs(layer) > 0
            nepochs = TrainEpochs(layer);
            train(layer, locs, nepochs, TrainMode); % locs is always 1, and TrainMode is only 2
        end
    end
    
    time2 = clock;
    disp(['Total Training Time: ', num2str(etime(time2, time1))]);
    
    if NET_SIZE == 32
        save('../VNanalysis/netTrained32.mat', 'NSYN', 'SQNET_SIZE', 'FiltImSize', 'SYNMATL1', 'SYNMATL2', 'SYNMATL3', 'SYNMATL4'); % this needs work
    elseif NET_SIZE == 64
        save('../VNanalysis/netTrained64.mat', 'NSYN', 'SQNET_SIZE', 'FiltImSize', 'SYNMATL1', 'SYNMATL2', 'SYNMATL3', 'SYNMATL4'); % this needs work
    elseif NET_SIZE == 128
        save('../VNanalysis/netTrained128.mat', 'NSYN', 'SQNET_SIZE', 'FiltImSize', 'SYNMATL1', 'SYNMATL2', 'SYNMATL3', 'SYNMATL4'); % this needs work
    elseif NET_SIZE == 256
        save('../VNanalysis/netTrained256.mat', 'NSYN', 'SQNET_SIZE', 'FiltImSize', 'SYNMATL1', 'SYNMATL2', 'SYNMATL3', 'SYNMATL4'); % this needs work
    end
    disp('Finished training');
end

if ischar(VNNetTrainedLoadName)
    disp(['Trained Net to be loaded is ' char(VNNetTrainedLoadName)]);
    load(char(VNNetTrainedLoadName));
end

% Include the following two lines only to test with interpolated view test images
% Interpolation = 1;  % Set this to 1 to test with interpolated images
% SetIMLIST(IMLIST, IMLIST_TEST, base_dir, Interpolation); % Set up the file names for the training and test objects in SetIMLIST.m

test(SweepFilename);    % present the stimuli, save the rates in EFRDATA, and save EFRDATA
info(NLAYERS);                       % single and multiple cell info analysis for layer 4
t2 = clock;
ElapsedTime = etime(t2, t1);
fprintf('Total runtime=%4.1f\n', ElapsedTime);

