Click here to Skip to main content
Click here to Skip to main content
Add your own
alternative version

Neural Network for Recognition of Handwritten Digits

, 5 Dec 2006
A convolutional neural network achieves 99.26% accuracy on a modified NIST database of hand-written digits.
// MNistDoc.cpp : implementation of the CMNistDoc class
//

#include "stdafx.h"
#include "MNist.h"

#include "MNistDoc.h"
#include "CntrItem.h"

extern CMNistApp theApp;

UINT CMNistDoc::m_iBackpropThreadIdentifier = 0;  // static member used by threads to identify themselves
UINT CMNistDoc::m_iTestingThreadIdentifier = 0;  

#include "SHLWAPI.H"	// for the path functions
#pragma comment( lib, "shlwapi.lib" )

#ifdef _DEBUG
#define new DEBUG_NEW
#undef THIS_FILE
static char THIS_FILE[] = __FILE__;
#endif

/////////////////////////////////////////////////////////////////////////////
// CMNistDoc

IMPLEMENT_DYNCREATE(CMNistDoc, COleDocument)

BEGIN_MESSAGE_MAP(CMNistDoc, COleDocument)
//{{AFX_MSG_MAP(CMNistDoc)
ON_BN_CLICKED(IDC_BUTTON_OPEN_MNIST_FILES, OnButtonOpenMnistFiles)
ON_BN_CLICKED(IDC_BUTTON_CLOSE_MNIST_FILES, OnButtonCloseMnistFiles)

//}}AFX_MSG_MAP
// Enable default OLE container implementation
ON_UPDATE_COMMAND_UI(ID_EDIT_PASTE, COleDocument::OnUpdatePasteMenu)
ON_UPDATE_COMMAND_UI(ID_EDIT_PASTE_LINK, COleDocument::OnUpdatePasteLinkMenu)
ON_UPDATE_COMMAND_UI(ID_OLE_EDIT_CONVERT, COleDocument::OnUpdateObjectVerbMenu)
ON_COMMAND(ID_OLE_EDIT_CONVERT, COleDocument::OnEditConvert)
ON_UPDATE_COMMAND_UI(ID_OLE_EDIT_LINKS, COleDocument::OnUpdateEditLinksMenu)
ON_COMMAND(ID_OLE_EDIT_LINKS, COleDocument::OnEditLinks)
ON_UPDATE_COMMAND_UI_RANGE(ID_OLE_VERB_FIRST, ID_OLE_VERB_LAST, COleDocument::OnUpdateObjectVerbMenu)
END_MESSAGE_MAP()

/////////////////////////////////////////////////////////////////////////////
// CMNistDoc construction/destruction

CMNistDoc::CMNistDoc()
{
	// Use OLE compound files
	EnableCompoundFile();
	
	// TODO: add one-time construction code here
	
	m_bFilesOpen = FALSE;
	m_bBackpropThreadAbortFlag = FALSE;
	m_bBackpropThreadsAreRunning = FALSE;
	m_cBackprops = 0;
	m_nAfterEveryNBackprops = 1;
	
	m_bTestingThreadsAreRunning = FALSE;
	m_bTestingThreadAbortFlag = FALSE;
	
	m_iNextTestingPattern = 0;
	m_iNextTrainingPattern = 0;
	
	::InitializeCriticalSection( &m_csTrainingPatterns );
	::InitializeCriticalSection( &m_csTestingPatterns );
	
	m_utxNeuralNet = ::CreateMutex( NULL, FALSE, NULL );  // anonymous mutex which is unowned initially
	
	ASSERT( m_utxNeuralNet != NULL );
	
	
	// allocate memory to store the distortion maps
	
	m_cCols = 29;
	m_cRows = 29;
	
	m_cCount = m_cCols * m_cRows;
	
	m_DispH = new double[ m_cCount ];
	m_DispV = new double[ m_cCount ];
	
	
	// create a gaussian kernel, which is constant, for use in generating elastic distortions
	
	int iiMid = GAUSSIAN_FIELD_SIZE/2;  // GAUSSIAN_FIELD_SIZE is strictly odd
	
	double twoSigmaSquared = 2.0 * (::GetPreferences().m_dElasticSigma) * (::GetPreferences().m_dElasticSigma);
	twoSigmaSquared = 1.0 /  twoSigmaSquared;
	double twoPiSigma = 1.0 / (::GetPreferences().m_dElasticSigma) * sqrt( 2.0 * 3.1415926535897932384626433832795 );
	
	for ( int col=0; col<GAUSSIAN_FIELD_SIZE; ++col )
	{
		for ( int row=0; row<GAUSSIAN_FIELD_SIZE; ++row )
		{
			m_GaussianKernel[ row ][ col ] = twoPiSigma * 
				( exp(- ( ((row-iiMid)*(row-iiMid) + (col-iiMid)*(col-iiMid)) * twoSigmaSquared ) ) );
		}
	}
	
	
	
}

CMNistDoc::~CMNistDoc()
{
	if ( m_bFilesOpen != FALSE )
	{
		CloseMnistFiles();
	}
	
	::DeleteCriticalSection( &m_csTrainingPatterns );
	::DeleteCriticalSection( &m_csTestingPatterns );
	
	::CloseHandle( m_utxNeuralNet );
	
	// delete memory of the distortion maps, allocated in constructor
	
	delete[] m_DispH;
	delete[] m_DispV;
	
	
}

//DEL BOOL CMNistDoc::OnOpenDocument(LPCTSTR lpszPathName) 
//DEL {
//DEL 	if (!COleDocument::OnOpenDocument(lpszPathName))
//DEL 		return FALSE;
//DEL 
//DEL 	
//DEL 	// TODO: Add your specialized creation code here
//DEL 	
//DEL 	return TRUE;
//DEL }

void CMNistDoc::DeleteContents() 
{
	// TODO: Add your specialized code here and/or call the base class
	
	COleDocument::DeleteContents();
	
	m_NN.Initialize();
	
}




BOOL CMNistDoc::OnNewDocument()
{
	if (!COleDocument::OnNewDocument())
		return FALSE;
	
	// TODO: add reinitialization code here
	// (SDI documents will reuse this document)
	
	// grab the mutex for the neural network
	
	CAutoMutex tlo( m_utxNeuralNet );
	
	// initialize and build the neural net
	
	NeuralNetwork& NN = m_NN;  // for easier nomenclature
	NN.Initialize();
	
	NNLayer* pLayer;
	
	int ii, jj, kk;
	int icNeurons = 0;
	int icWeights = 0;
	double initWeight;
	CString label;
	
	// layer zero, the input layer.
	// Create neurons: exactly the same number of neurons as the input
	// vector of 29x29=841 pixels, and no weights/connections
	
	pLayer = new NNLayer( _T("Layer00") );
	NN.m_Layers.push_back( pLayer );
	
	for ( ii=0; ii<841; ++ii )
	{
		label.Format( _T("Layer00_Neuron%04d_Num%06d"), ii, icNeurons );
		pLayer->m_Neurons.push_back( new NNNeuron( (LPCTSTR)label ) );
		icNeurons++;
	}
	
#define UNIFORM_PLUS_MINUS_ONE ( (double)(2.0 * rand())/RAND_MAX - 1.0 )
	
	// layer one:
	// This layer is a convolutional layer that has 6 feature maps.  Each feature 
	// map is 13x13, and each unit in the feature maps is a 5x5 convolutional kernel
	// of the input layer.
	// So, there are 13x13x6 = 1014 neurons, (5x5+1)x6 = 156 weights
	
	pLayer = new NNLayer( _T("Layer01"), pLayer );
	NN.m_Layers.push_back( pLayer );
	
	for ( ii=0; ii<1014; ++ii )
	{
		label.Format( _T("Layer01_Neuron%04d_Num%06d"), ii, icNeurons );
		pLayer->m_Neurons.push_back( new NNNeuron( (LPCTSTR)label ) );
		icNeurons++;
	}
	
	for ( ii=0; ii<156; ++ii )
	{
		label.Format( _T("Layer01_Weight%04d_Num%06d"), ii, icWeights );
		initWeight = 0.05 * UNIFORM_PLUS_MINUS_ONE;
		pLayer->m_Weights.push_back( new NNWeight( (LPCTSTR)label, initWeight ) );
	}
	
	// interconnections with previous layer: this is difficult
	// The previous layer is a top-down bitmap image that has been padded to size 29x29
	// Each neuron in this layer is connected to a 5x5 kernel in its feature map, which 
	// is also a top-down bitmap of size 13x13.  We move the kernel by TWO pixels, i.e., we
	// skip every other pixel in the input image
	
	int kernelTemplate[25] = {
		0,  1,  2,  3,  4,
		29, 30, 31, 32, 33,
		58, 59, 60, 61, 62,
		87, 88, 89, 90, 91,
		116,117,118,119,120 };
		
	int iNumWeight;
		
	int fm;
		
	for ( fm=0; fm<6; ++fm)
	{
		for ( ii=0; ii<13; ++ii )
		{
			for ( jj=0; jj<13; ++jj )
			{
				iNumWeight = fm * 26;  // 26 is the number of weights per feature map
				NNNeuron& n = *( pLayer->m_Neurons[ jj + ii*13 + fm*169 ] );
				
				n.AddConnection( ULONG_MAX, iNumWeight++ );  // bias weight
				
				for ( kk=0; kk<25; ++kk )
				{
					// note: max val of index == 840, corresponding to 841 neurons in prev layer
					n.AddConnection( 2*jj + 58*ii + kernelTemplate[kk], iNumWeight++ );
				}
			}
		}
	}
	
	
	// layer two:
	// This layer is a convolutional layer that has 50 feature maps.  Each feature 
	// map is 5x5, and each unit in the feature maps is a 5x5 convolutional kernel
	// of corresponding areas of all 6 of the previous layers, each of which is a 13x13 feature map
	// So, there are 5x5x50 = 1250 neurons, (5x5+1)x6x50 = 7800 weights
	
	pLayer = new NNLayer( _T("Layer02"), pLayer );
	NN.m_Layers.push_back( pLayer );
	
	for ( ii=0; ii<1250; ++ii )
	{
		label.Format( _T("Layer02_Neuron%04d_Num%06d"), ii, icNeurons );
		pLayer->m_Neurons.push_back( new NNNeuron( (LPCTSTR)label ) );
		icNeurons++;
	}
	
	for ( ii=0; ii<7800; ++ii )
	{
		label.Format( _T("Layer02_Weight%04d_Num%06d"), ii, icWeights );
		initWeight = 0.05 * UNIFORM_PLUS_MINUS_ONE;
		pLayer->m_Weights.push_back( new NNWeight( (LPCTSTR)label, initWeight ) );
	}
	
	// Interconnections with previous layer: this is difficult
	// Each feature map in the previous layer is a top-down bitmap image whose size
	// is 13x13, and there are 6 such feature maps.  Each neuron in one 5x5 feature map of this 
	// layer is connected to a 5x5 kernel positioned correspondingly in all 6 parent
	// feature maps, and there are individual weights for the six different 5x5 kernels.  As
	// before, we move the kernel by TWO pixels, i.e., we
	// skip every other pixel in the input image.  The result is 50 different 5x5 top-down bitmap
	// feature maps
	
	int kernelTemplate2[25] = {
		0,  1,  2,  3,  4,
		13, 14, 15, 16, 17, 
		26, 27, 28, 29, 30,
		39, 40, 41, 42, 43, 
		52, 53, 54, 55, 56   };
		
		
	for ( fm=0; fm<50; ++fm)
	{
		for ( ii=0; ii<5; ++ii )
		{
			for ( jj=0; jj<5; ++jj )
			{
				iNumWeight = fm * 26;  // 26 is the number of weights per feature map
				NNNeuron& n = *( pLayer->m_Neurons[ jj + ii*5 + fm*25 ] );
				
				n.AddConnection( ULONG_MAX, iNumWeight++ );  // bias weight
				
				for ( kk=0; kk<25; ++kk )
				{
					// note: max val of index == 1013, corresponding to 1014 neurons in prev layer
					n.AddConnection(       2*jj + 26*ii + kernelTemplate2[kk], iNumWeight++ );
					n.AddConnection( 169 + 2*jj + 26*ii + kernelTemplate2[kk], iNumWeight++ );
					n.AddConnection( 338 + 2*jj + 26*ii + kernelTemplate2[kk], iNumWeight++ );
					n.AddConnection( 507 + 2*jj + 26*ii + kernelTemplate2[kk], iNumWeight++ );
					n.AddConnection( 676 + 2*jj + 26*ii + kernelTemplate2[kk], iNumWeight++ );
					n.AddConnection( 845 + 2*jj + 26*ii + kernelTemplate2[kk], iNumWeight++ );
				}
			}
		}
	}
			
	
	// layer three:
	// This layer is a fully-connected layer with 100 units.  Since it is fully-connected,
	// each of the 100 neurons in the layer is connected to all 1250 neurons in
	// the previous layer.
	// So, there are 100 neurons and 100*(1250+1)=125100 weights
	
	pLayer = new NNLayer( _T("Layer03"), pLayer );
	NN.m_Layers.push_back( pLayer );
	
	for ( ii=0; ii<100; ++ii )
	{
		label.Format( _T("Layer03_Neuron%04d_Num%06d"), ii, icNeurons );
		pLayer->m_Neurons.push_back( new NNNeuron( (LPCTSTR)label ) );
		icNeurons++;
	}
	
	for ( ii=0; ii<125100; ++ii )
	{
		label.Format( _T("Layer03_Weight%04d_Num%06d"), ii, icWeights );
		initWeight = 0.05 * UNIFORM_PLUS_MINUS_ONE;
		pLayer->m_Weights.push_back( new NNWeight( (LPCTSTR)label, initWeight ) );
	}
	
	// Interconnections with previous layer: fully-connected
	
	iNumWeight = 0;  // weights are not shared in this layer
	
	for ( fm=0; fm<100; ++fm )
	{
		NNNeuron& n = *( pLayer->m_Neurons[ fm ] );
		n.AddConnection( ULONG_MAX, iNumWeight++ );  // bias weight
		
		for ( ii=0; ii<1250; ++ii )
		{
			n.AddConnection( ii, iNumWeight++ );
		}
	}
	
			
			
	// layer four, the final (output) layer:
	// This layer is a fully-connected layer with 10 units.  Since it is fully-connected,
	// each of the 10 neurons in the layer is connected to all 100 neurons in
	// the previous layer.
	// So, there are 10 neurons and 10*(100+1)=1010 weights
	
	pLayer = new NNLayer( _T("Layer04"), pLayer );
	NN.m_Layers.push_back( pLayer );
	
	for ( ii=0; ii<10; ++ii )
	{
		label.Format( _T("Layer04_Neuron%04d_Num%06d"), ii, icNeurons );
		pLayer->m_Neurons.push_back( new NNNeuron( (LPCTSTR)label ) );
		icNeurons++;
	}
	
	for ( ii=0; ii<1010; ++ii )
	{
		label.Format( _T("Layer04_Weight%04d_Num%06d"), ii, icWeights );
		initWeight = 0.05 * UNIFORM_PLUS_MINUS_ONE;
		pLayer->m_Weights.push_back( new NNWeight( (LPCTSTR)label, initWeight ) );
	}
	
	// Interconnections with previous layer: fully-connected
	
	iNumWeight = 0;  // weights are not shared in this layer
	
	for ( fm=0; fm<10; ++fm )
	{
		NNNeuron& n = *( pLayer->m_Neurons[ fm ] );
		n.AddConnection( ULONG_MAX, iNumWeight++ );  // bias weight
		
		for ( ii=0; ii<100; ++ii )
		{
			n.AddConnection( ii, iNumWeight++ );
		}
	}
	
	
	SetModifiedFlag( TRUE );
	
	return TRUE;
}



/////////////////////////////////////////////////////////////////////////////
// CMNistDoc serialization

void CMNistDoc::Serialize(CArchive& ar)
{
	if (ar.IsStoring())
	{
		// TODO: add storing code here
		
	}
	else
	{
		// TODO: add loading code here
		
	}
	
	{
		// grab the mutex for the neural network in a local scope, and then serialize
		
		CAutoMutex tlo( m_utxNeuralNet );
		
		m_NN.Serialize( ar );
		
	}
	
	// Calling the base class COleDocument enables serialization
	//  of the container document's COleClientItem objects.
	COleDocument::Serialize(ar);
}

/////////////////////////////////////////////////////////////////////////////
// CMNistDoc diagnostics

#ifdef _DEBUG
void CMNistDoc::AssertValid() const
{
	COleDocument::AssertValid();
}

void CMNistDoc::Dump(CDumpContext& dc) const
{
	COleDocument::Dump(dc);
}
#endif //_DEBUG

/////////////////////////////////////////////////////////////////////////////
// CMNistDoc commands

void CMNistDoc::OnButtonOpenMnistFiles() 
{
	// opens the MNIST image data file and the label file
	
	struct FILEBEGINNING
	{
		union
		{
			char raw0[ 4 ];
			int nMagic;
		};
		union
		{
			char raw1[ 4 ];
			UINT nItems;
		};
		union
		{
			char raw2[ 4 ];
			int nRows;
		};
		union
		{
			char raw3[ 4 ];
			int nCols;
		};
	};
	
	if ( m_bFilesOpen != FALSE )
	{
		::MessageBox( NULL, _T("Files already open"), _T("Files already open"), MB_ICONEXCLAMATION );
		return;
	}
	
	CWaitCursor wc;
	
	// we need to open four files: (1) training images, (2) training labels, (3) testing images, (4) testing labels
	
	
	// (1) Training images
	
	CFileDialog fd( TRUE );  // constructs an "open" dialog; we are happy with the remaining defaults
	
	fd.m_ofn.lpstrFilter = _T("All Files (*.*)\0*.*\0\0");
	fd.m_ofn.lpstrTitle = _T("TRAINING Images");
	fd.m_ofn.lpstrInitialDir = theApp.m_sModulePath;
	
	if ( fd.DoModal() != IDOK ) return;
	
	// store chosen directory so that subsequent open-file dialogs can use it
	
	CString newPath = fd.GetPathName();
	::PathMakePretty( newPath.GetBuffer(255) );
	::PathRemoveFileSpec( newPath.GetBuffer(255) );
	newPath.ReleaseBuffer();
	
	// open the file
	
	CFileException fe;
	
	if ( m_fileTrainingImages.Open( (LPCTSTR)fd.GetPathName(), CFile::modeRead | CFile::shareDenyNone, &fe ) != 0 )
	{
		// opened successfully
		// check for expected magic number and for expected dimensionality
		
		FILEBEGINNING fileBegin;
		m_fileTrainingImages.Read( fileBegin.raw0, sizeof( fileBegin ) );
		
		// convert endian-ness (using winsock functions for convenience)
		
		fileBegin.nMagic = ::ntohl( fileBegin.nMagic );
		fileBegin.nItems = ::ntohl( fileBegin.nItems );
		fileBegin.nRows = ::ntohl( fileBegin.nRows );
		fileBegin.nCols = ::ntohl( fileBegin.nCols );
		
		// check against expected values
		
		if ( fileBegin.nMagic != ::GetPreferences().m_nMagicTrainingImages ||
			fileBegin.nItems != ::GetPreferences().m_nItemsTrainingImages ||
			fileBegin.nRows != ::GetPreferences().m_nRowsImages ||
			fileBegin.nCols != ::GetPreferences().m_nColsImages )
		{
			// file is not configured as expected
			CString msg;
			msg.Format( _T("File header for training images contains unexpected values as follows\n")
				_T("Magic number: got %i, expected %i \n")
				_T("Number of items: got %i, expected %i \n")
				_T("Number of rows: got %i, expected %i \n")
				_T("Number of columns: got %i, expected %i \n")
				_T("\nCheck integrity of file for training images and/or adjust the expected values in the INI file"),
				fileBegin.nMagic, ::GetPreferences().m_nMagicTrainingImages,
				fileBegin.nItems, ::GetPreferences().m_nItemsTrainingImages,
				fileBegin.nRows, ::GetPreferences().m_nRowsImages,
				fileBegin.nCols, ::GetPreferences().m_nColsImages );
			
			::MessageBox( NULL, msg, _T("Unexpected Format"), MB_OK|MB_ICONEXCLAMATION );
			
			// close all files
			
			m_fileTrainingImages.Close();
			
			return;
		}
	}
	else
	{
		// could not open image file
		
		TRACE(_T("%s(%i): Could not open file for training images \n"), __FILE__,__LINE__);
		
		CString msg;
		TCHAR szCause[255] = {0};
		
		fe.GetErrorMessage(szCause, 254);
		
		msg.Format( _T("Image file for training images could not be opened\n")
			_T("File name: %s\nReason: %s\nError code: %d"), 
			fe.m_strFileName, szCause, fe.m_cause );
		
		::MessageBox( NULL, msg, _T("Failed to open file for training images"), MB_OK|MB_ICONEXCLAMATION );
		
		return;
	}
	
	
	
	// (2) Training labels
	
	CFileDialog fd2( TRUE );  // constructs an "open" dialog; we are happy with the remaining defaults
	
	fd2.m_ofn.lpstrFilter = _T("All Files (*.*)\0*.*\0\0");
	fd2.m_ofn.lpstrTitle = _T("TRAINING Labels");
	fd2.m_ofn.lpstrInitialDir = newPath;
	
	if ( fd2.DoModal() != IDOK ) 
	{
		// close the images file too
		
		m_fileTrainingImages.Close();
		return;
	}
	
	
	// open the file
	
	if ( m_fileTrainingLabels.Open( (LPCTSTR)fd2.GetPathName(), CFile::modeRead | CFile::shareDenyNone, &fe ) != 0 )
	{
		// opened successfully
		// check for expected magic number and for expected dimensionality
		
		FILEBEGINNING fileBegin;
		m_fileTrainingLabels.Read( fileBegin.raw0, sizeof( fileBegin ) );
		
		// convert endian-ness (using winsock functions for convenience)
		
		fileBegin.nMagic = ::ntohl( fileBegin.nMagic );
		fileBegin.nItems = ::ntohl( fileBegin.nItems );
		
		// check against expected values
		
		if ( fileBegin.nMagic != ::GetPreferences().m_nMagicTrainingLabels ||
			fileBegin.nItems != ::GetPreferences().m_nItemsTrainingLabels )
		{
			// file is not configured as expected
			CString msg;
			msg.Format( _T("File header for training labels contains unexpected values as follows\n")
				_T("Magic number: got %i, expected %i \n")
				_T("Number of items: got %i, expected %i \n")
				_T("\nCheck integrity of file for training labels and/or adjust the expected values in the INI file"),
				fileBegin.nMagic, ::GetPreferences().m_nMagicTrainingLabels,
				fileBegin.nItems, ::GetPreferences().m_nItemsTrainingLabels );
			
			::MessageBox( NULL, msg, _T("Unexpected Format"), MB_OK|MB_ICONEXCLAMATION );
			
			// close all files
			
			m_fileTrainingImages.Close();
			m_fileTrainingLabels.Close();
			
			return;
		}
	}
	else
	{
		// could not open label file
		
		TRACE(_T("%s(%i): Could not open file for training labels \n"), __FILE__,__LINE__);
		
		CString msg;
		TCHAR szCause[255] = {0};
		
		fe.GetErrorMessage(szCause, 254);
		
		msg.Format( _T("Label file for training labels could not be opened\n")
			_T("File name: %s\nReason: %s\nError code: %d"), 
			fe.m_strFileName, szCause, fe.m_cause );
		
		::MessageBox( NULL, msg, _T("Failed to open file for training labels"), MB_OK|MB_ICONEXCLAMATION );
		
		
		// close the already-opened files too
		
		m_fileTrainingImages.Close();
		
		return;
	}
	
	
	
	// (3) Testing images
	
	CFileDialog fd3( TRUE );  // constructs an "open" dialog; we are happy with the remaining defaults
	
	fd3.m_ofn.lpstrFilter = _T("All Files (*.*)\0*.*\0\0");
	fd3.m_ofn.lpstrTitle = _T("TESTING Images");
	fd3.m_ofn.lpstrInitialDir = newPath;
	
	if ( fd3.DoModal() != IDOK ) 
	{
		// close the already-opened files too
		
		m_fileTrainingLabels.Close();
		m_fileTrainingImages.Close();
		
		return;
	}
	
	
	// open the file
	
	if ( m_fileTestingImages.Open( (LPCTSTR)fd3.GetPathName(), CFile::modeRead | CFile::shareDenyNone, &fe ) != 0 )
	{
		// opened successfully
		// check for expected magic number and for expected dimensionality
		
		FILEBEGINNING fileBegin;
		m_fileTestingImages.Read( fileBegin.raw0, sizeof( fileBegin ) );
		
		// convert endian-ness (using winsock functions for convenience)
		
		fileBegin.nMagic = ::ntohl( fileBegin.nMagic );
		fileBegin.nItems = ::ntohl( fileBegin.nItems );
		fileBegin.nRows = ::ntohl( fileBegin.nRows );
		fileBegin.nCols = ::ntohl( fileBegin.nCols );
		
		// check against expected values
		
		if ( fileBegin.nMagic != ::GetPreferences().m_nMagicTestingImages ||
			fileBegin.nItems != ::GetPreferences().m_nItemsTestingImages ||
			fileBegin.nRows != ::GetPreferences().m_nRowsImages ||
			fileBegin.nCols != ::GetPreferences().m_nColsImages )
		{
			// file is not configured as expected
			CString msg;
			msg.Format( _T("File header for testing images contains unexpected values as follows\n")
				_T("Magic number: got %i, expected %i \n")
				_T("Number of items: got %i, expected %i \n")
				_T("Number of rows: got %i, expected %i \n")
				_T("Number of columns: got %i, expected %i \n")
				_T("\nCheck integrity of file for testing images and/or adjust the expected values in the INI file"),
				fileBegin.nMagic, ::GetPreferences().m_nMagicTestingImages,
				fileBegin.nItems, ::GetPreferences().m_nItemsTestingImages,
				fileBegin.nRows, ::GetPreferences().m_nRowsImages,
				fileBegin.nCols, ::GetPreferences().m_nColsImages );
			
			::MessageBox( NULL, msg, _T("Unexpected Format"), MB_OK|MB_ICONEXCLAMATION );
			
			
			// close all files
			
			m_fileTestingImages.Close();
			m_fileTrainingLabels.Close();
			m_fileTrainingImages.Close();
			
			return;
		}
	}
	else
	{
		// could not open image file
		
		TRACE(_T("%s(%i): Could not open file for testing images \n"), __FILE__,__LINE__);
		
		CString msg;
		TCHAR szCause[255] = {0};
		
		fe.GetErrorMessage(szCause, 254);
		
		msg.Format( _T("Image file for testing images could not be opened\n")
			_T("File name: %s\nReason: %s\nError code: %d"), 
			fe.m_strFileName, szCause, fe.m_cause );
		
		::MessageBox( NULL, msg, _T("Failed to open file for testing images"), MB_OK|MB_ICONEXCLAMATION );
		
		
		// close the already-opened files too
		
		m_fileTrainingLabels.Close();
		m_fileTrainingImages.Close();
		
		return;
	}
	
	
	
	// (4) Testing labels
	
	CFileDialog fd4( TRUE );  // constructs an "open" dialog; we are happy with the remaining defaults
	
	fd4.m_ofn.lpstrFilter = _T("All Files (*.*)\0*.*\0\0");
	fd4.m_ofn.lpstrTitle = _T("TESTING Labels");
	fd4.m_ofn.lpstrInitialDir = newPath;
	
	if ( fd4.DoModal() != IDOK ) 
	{
		// close the already-opened files too
		
		m_fileTestingImages.Close();
		m_fileTrainingLabels.Close();
		m_fileTrainingImages.Close();
		
		return;
	}
	
	
	// open the file
	
	if ( m_fileTestingLabels.Open( (LPCTSTR)fd4.GetPathName(), CFile::modeRead | CFile::shareDenyNone, &fe ) != 0 )
	{
		// opened successfully
		// check for expected magic number and for expected dimensionality
		
		FILEBEGINNING fileBegin;
		m_fileTestingLabels.Read( fileBegin.raw0, sizeof( fileBegin ) );
		
		// convert endian-ness (using winsock functions for convenience)
		
		fileBegin.nMagic = ::ntohl( fileBegin.nMagic );
		fileBegin.nItems = ::ntohl( fileBegin.nItems );
		
		// check against expected values
		
		if ( fileBegin.nMagic != ::GetPreferences().m_nMagicTestingLabels ||
			fileBegin.nItems != ::GetPreferences().m_nItemsTestingLabels )
		{
			// file is not configured as expected
			CString msg;
			msg.Format( _T("File header for testing labels contains unexpected values as follows\n")
				_T("Magic number: got %i, expected %i \n")
				_T("Number of items: got %i, expected %i \n")
				_T("\nCheck integrity of file for testing labels and/or adjust the expected values in the INI file"),
				fileBegin.nMagic, ::GetPreferences().m_nMagicTestingLabels,
				fileBegin.nItems, ::GetPreferences().m_nItemsTestingLabels );
			
			::MessageBox( NULL, msg, _T("Unexpected Format"), MB_OK|MB_ICONEXCLAMATION );
			
			// close all files
			
			m_fileTestingLabels.Close();
			m_fileTestingImages.Close();
			m_fileTrainingImages.Close();
			m_fileTrainingLabels.Close();
			
			return;
		}
	}
	else
	{
		// could not open label file
		
		TRACE(_T("%s(%i): Could not open file for testing labels \n"), __FILE__,__LINE__);
		
		CString msg;
		TCHAR szCause[255] = {0};
		
		fe.GetErrorMessage(szCause, 254);
		
		msg.Format( _T("Label file for testing labels could not be opened\n")
			_T("File name: %s\nReason: %s\nError code: %d"), 
			fe.m_strFileName, szCause, fe.m_cause );
		
		::MessageBox( NULL, msg, _T("Failed to open file for testing labels"), MB_OK|MB_ICONEXCLAMATION );
		
		// close all already-opened files
		
		m_fileTestingImages.Close();
		m_fileTrainingImages.Close();
		m_fileTrainingLabels.Close();
		
		return;
	}
	
	
	
	// all four files are opened, and all four contain expected header information
	
	m_bFilesOpen = TRUE;
	
	ASSERT( g_cImageSize == 28 );
	
}

void CMNistDoc::OnButtonCloseMnistFiles()
{
	CloseMnistFiles();
}

void CMNistDoc::CloseMnistFiles()
{
	m_fileTestingImages.Close();
	m_fileTestingLabels.Close();
	m_fileTrainingImages.Close();
	m_fileTrainingLabels.Close();
	
	m_bFilesOpen = FALSE;
	
}



double CMNistDoc::GetCurrentEta()
{
	return m_NN.m_etaLearningRate;
}


double CMNistDoc::GetPreviousEta()
{
	// provided because threads might change the current eta before we are able to read it
	
	return m_NN.m_etaLearningRatePrevious;
}


UINT CMNistDoc::GetCurrentTrainingPatternNumber( BOOL bFromRandomizedPatternSequence /* =FALSE */ )
{
	// returns the current number of the training pattern, either from the straight sequence, or from
	// the randomized sequence
	
	UINT iRet;
	
	if ( bFromRandomizedPatternSequence == FALSE )
	{
		iRet = m_iNextTrainingPattern;
	}
	else
	{
		iRet = m_iRandomizedTrainingPatternSequence[ m_iNextTrainingPattern ];
	}
	
	return iRet;
}



// UNIFORM_ZERO_THRU_ONE gives a uniformly-distributed number between zero (inclusive) and one (exclusive)

#define UNIFORM_ZERO_THRU_ONE ( (double)(rand())/(RAND_MAX + 1 ) ) 


void CMNistDoc::RandomizeTrainingPatternSequence()
{
	// randomizes the order of m_iRandomizedTrainingPatternSequence, which is a UINT array
	// holding the numbers 0..59999 in random order
	
	CAutoCS tlo( m_csTrainingPatterns );
	
	UINT ii, jj, iiMax, iiTemp;
	
	iiMax = ::GetPreferences().m_nItemsTrainingImages;
	
	ASSERT( iiMax == 60000 );  // requirement of sloppy and unimaginative code
	
	// initialize array in sequential order
	
	for ( ii=0; ii<iiMax; ++ii )
	{
		m_iRandomizedTrainingPatternSequence[ ii ] = ii;  
	}
	
	
	// now at each position, swap with a random position
	
	for ( ii=0; ii<iiMax; ++ii )
	{
		jj = (UINT)( UNIFORM_ZERO_THRU_ONE * iiMax );
		
		ASSERT( jj < iiMax );
		
		iiTemp = m_iRandomizedTrainingPatternSequence[ ii ];
		m_iRandomizedTrainingPatternSequence[ ii ] = m_iRandomizedTrainingPatternSequence[ jj ];
		m_iRandomizedTrainingPatternSequence[ jj ] = iiTemp;
	}
	
}



UINT CMNistDoc::GetNextTrainingPattern(unsigned char *pArray /* =NULL */, int *pLabel /* =NULL */, 
									   BOOL bFlipGrayscale /* =TRUE */, BOOL bFromRandomizedPatternSequence /* =TRUE */, UINT* iSequenceNum /* =NULL */)
{
	// returns the number of the pattern corresponding to the pattern that will be stored in pArray
	// if BOOL bFromRandomizedPatternSequence is TRUE (which is the default) then the pattern
	// stored will be a pattern from the randomized sequence; otherwise the pattern will be a straight
	// sequential run through all the training patterns, from start to finish.  The sequence number,
	// which runs from 0..59999 monotonically, is returned in iSequenceNum (if it's not NULL)
	
	CAutoCS tlo( m_csTrainingPatterns );
	
	UINT iPatternNum;
	
	if ( bFromRandomizedPatternSequence == FALSE )
	{
		iPatternNum = m_iNextTrainingPattern;
	}
	else
	{
		iPatternNum = m_iRandomizedTrainingPatternSequence[ m_iNextTrainingPattern ];
	}

	ASSERT( iPatternNum < ::GetPreferences().m_nItemsTrainingImages );
	
	GetTrainingPatternArrayValues( iPatternNum, pArray, pLabel, bFlipGrayscale );
	
	if ( iSequenceNum != NULL )
	{
		*iSequenceNum = m_iNextTrainingPattern;
	}
	
	m_iNextTrainingPattern++;
	
	if ( m_iNextTrainingPattern >= ::GetPreferences().m_nItemsTrainingImages )
	{
		m_iNextTrainingPattern = 0;
	}
	
	return iPatternNum;
}




UINT CMNistDoc::GetRandomTrainingPattern(unsigned char *pArray /* =NULL */, int *pLabel /* =NULL */, BOOL bFlipGrayscale /* =TRUE */ )
{
	// returns the number of the pattern corresponding to the pattern stored in pArray
	
	CAutoCS tlo( m_csTrainingPatterns );
	
	UINT patternNum = (UINT)( UNIFORM_ZERO_THRU_ONE * (::GetPreferences().m_nItemsTrainingImages - 1) );
	
	GetTrainingPatternArrayValues( patternNum, pArray, pLabel, bFlipGrayscale );
	
	return patternNum;
}





void CMNistDoc::GetTrainingPatternArrayValues(int iNumImage /* =0 */, unsigned char *pArray /* =NULL */, int *pLabel /* =NULL */,
											  BOOL bFlipGrayscale /* =TRUE */ )
{
	// fills an unsigned char array with gray values, corresponding to iNumImage, and also
	// returns the label for the image
	
	CAutoCS tlo( m_csTrainingPatterns );
	
	int cCount = g_cImageSize*g_cImageSize;
	int fPos;
	
	if ( m_bFilesOpen != FALSE )
	{
		if ( pArray != NULL )
		{
			fPos = 16 + iNumImage*cCount;  // 16 compensates for file header info
			m_fileTrainingImages.Seek( fPos, CFile::begin );
			m_fileTrainingImages.Read( pArray, cCount );
			
			if ( bFlipGrayscale != FALSE )
			{
				for ( int ii=0; ii<cCount; ++ii )
				{
					pArray[ ii ] = 255 - pArray[ ii ];
				}
			}
		}
		
		if ( pLabel != NULL )
		{
			fPos = 8 + iNumImage;
			char r;
			m_fileTrainingLabels.Seek( fPos, CFile::begin );
			m_fileTrainingLabels.Read( &r, 1 );  // single byte
			
			*pLabel = r;
		}
	}
	else  // no files are open: return a simple gray wedge
	{
		if ( pArray != NULL )
		{
			for ( int ii=0; ii<cCount; ++ii )
			{
				pArray[ ii ] = ii*255/cCount;
			}
		}
		
		if ( pLabel != NULL )
		{
			*pLabel = INT_MAX;
		}
	}
}



UINT CMNistDoc::GetNextTestingPatternNumber()
{
	return m_iNextTestingPattern;
}



UINT CMNistDoc::GetNextTestingPattern(unsigned char *pArray /* =NULL */, int *pLabel /* =NULL */, BOOL bFlipGrayscale /* =TRUE */ )
{
	// returns the number of the pattern corresponding to the pattern stored in pArray
	
	CAutoCS tlo( m_csTestingPatterns );
	
	
	GetTestingPatternArrayValues( m_iNextTestingPattern, pArray, pLabel, bFlipGrayscale );
	
	UINT iRet = m_iNextTestingPattern;
	m_iNextTestingPattern++;
	
	if ( m_iNextTestingPattern >= ::GetPreferences().m_nItemsTestingImages )
	{
		m_iNextTestingPattern = 0;
	}
	
	return iRet ;
}


void CMNistDoc::GetTestingPatternArrayValues(int iNumImage /* =0 */, unsigned char *pArray /* =NULL */, int *pLabel /* =NULL */,
											 BOOL bFlipGrayscale /* =TRUE */ )
{
	// fills an unsigned char array with gray values, corresponding to iNumImage, and also
	// returns the label for the image
	
	CAutoCS tlo( m_csTestingPatterns );
	
	int cCount = g_cImageSize*g_cImageSize;
	int fPos;
	
	if ( m_bFilesOpen != FALSE )
	{
		if ( pArray != NULL )
		{
			fPos = 16 + iNumImage*cCount;  // 16 compensates for file header info
			m_fileTestingImages.Seek( fPos, CFile::begin );
			m_fileTestingImages.Read( pArray, cCount );
			
			if ( bFlipGrayscale != FALSE )
			{
				for ( int ii=0; ii<cCount; ++ii )
				{
					pArray[ ii ] = 255 - pArray[ ii ];
				}
			}
		}
		
		if ( pLabel != NULL )
		{
			fPos = 8 + iNumImage;
			char r;
			m_fileTestingLabels.Seek( fPos, CFile::begin );
			m_fileTestingLabels.Read( &r, 1 );  // single byte
			
			*pLabel = r;
		}
	}
	else  // no files are open: return a simple gray wedge
	{
		if ( pArray != NULL )
		{
			for ( int ii=0; ii<cCount; ++ii )
			{
				pArray[ ii ] = ii*255/cCount;
			}
		}
		
		if ( pLabel != NULL )
		{
			*pLabel = INT_MAX;
		}
	}
}



void CMNistDoc::GenerateDistortionMap( double severityFactor /* =1.0 */ )
{
	// generates distortion maps in each of the horizontal and vertical directions
	// Three distortions are applied: a scaling, a rotation, and an elastic distortion
	// Since these are all linear tranformations, we can simply add them together, after calculation
	// one at a time
	
	// The input parameter, severityFactor, let's us control the severity of the distortions relative
	// to the default values.  For example, if we only want half as harsh a distortion, set
	// severityFactor == 0.5
	
	// First, elastic distortion, per Patrice Simard, "Best Practices For Convolutional Neural Networks..."
	// at page 2.
	// Three-step process: seed array with uniform randoms, filter with a gaussian kernel, normalize (scale)
	
	int row, col;
	double* uniformH = new double[ m_cCount ];
	double* uniformV = new double[ m_cCount ];
	
	
	for ( col=0; col<m_cCols; ++col )
	{
		for ( row=0; row<m_cRows; ++row )
		{
			At( uniformH, row, col ) = UNIFORM_PLUS_MINUS_ONE;
			At( uniformV, row, col ) = UNIFORM_PLUS_MINUS_ONE;
		}
	}
	
	// filter with gaussian
	
	double fConvolvedH, fConvolvedV;
	double fSampleH, fSampleV;
	double elasticScale = severityFactor * ::GetPreferences().m_dElasticScaling;
	int xxx, yyy, xxxDisp, yyyDisp;
	int iiMid = GAUSSIAN_FIELD_SIZE/2;  // GAUSSIAN_FIELD_SIZE is strictly odd
	
	for ( col=0; col<m_cCols; ++col )
	{
		for ( row=0; row<m_cRows; ++row )
		{
			fConvolvedH = 0.0;
			fConvolvedV = 0.0;
			
			for ( xxx=0; xxx<GAUSSIAN_FIELD_SIZE; ++xxx )
			{
				for ( yyy=0; yyy<GAUSSIAN_FIELD_SIZE; ++yyy )
				{
					xxxDisp = col - iiMid + xxx;
					yyyDisp = row - iiMid + yyy;
					
					if ( xxxDisp<0 || xxxDisp>=m_cCols || yyyDisp<0 || yyyDisp>=m_cRows )
					{
						fSampleH = 0.0;
						fSampleV = 0.0;
					}
					else
					{
						fSampleH = At( uniformH, yyyDisp, xxxDisp );
						fSampleV = At( uniformV, yyyDisp, xxxDisp );
					}
					
					fConvolvedH += fSampleH * m_GaussianKernel[ yyy ][ xxx ];
					fConvolvedV += fSampleV * m_GaussianKernel[ yyy ][ xxx ];
				}
			}
			
			At( m_DispH, row, col ) = elasticScale * fConvolvedH;
			At( m_DispV, row, col ) = elasticScale * fConvolvedV;
		}
	}
	
	delete[] uniformH;
	delete[] uniformV;
	
	// next, the scaling of the image by a random scale factor
	// Horizontal and vertical directions are scaled independently
	
	double dSFHoriz = severityFactor * ::GetPreferences().m_dMaxScaling / 100.0 * UNIFORM_PLUS_MINUS_ONE;  // m_dMaxScaling is a percentage
	double dSFVert = severityFactor * ::GetPreferences().m_dMaxScaling / 100.0 * UNIFORM_PLUS_MINUS_ONE;  // m_dMaxScaling is a percentage

	
	int iMid = m_cRows/2;
	
	for ( row=0; row<m_cRows; ++row )
	{
		for ( col=0; col<m_cCols; ++col )
		{
			At( m_DispH, row, col ) += dSFHoriz * ( col-iMid );
			At( m_DispV, row, col ) -= dSFVert * ( iMid-row );  // negative because of top-down bitmap
		}
	}
	
	
	// finally, apply a rotation
	
	double angle = severityFactor * ::GetPreferences().m_dMaxRotation * UNIFORM_PLUS_MINUS_ONE;
	angle = angle * 3.1415926535897932384626433832795 / 180.0;  // convert from degrees to radians
	
	double cosAngle = cos( angle );
	double sinAngle = sin( angle );
	
	for ( row=0; row<m_cRows; ++row )
	{
		for ( col=0; col<m_cCols; ++col )
		{
			At( m_DispH, row, col ) += ( col-iMid ) * ( cosAngle - 1 ) - ( iMid-row ) * sinAngle;
			At( m_DispV, row, col ) -= ( iMid-row ) * ( cosAngle - 1 ) + ( col-iMid ) * sinAngle;  // negative because of top-down bitmap
		}
	}
	
}


void CMNistDoc::ApplyDistortionMap(double *inputVector)
{
	// applies the current distortion map to the input vector
	
	// For the mapped array, we assume that 0.0 == background, and 1.0 == full intensity information
	// This is different from the input vector, in which +1.0 == background (white), and 
	// -1.0 == information (black), so we must convert one to the other
	
	std::vector< std::vector< double > >   mappedVector( m_cRows, std::vector< double >( m_cCols, 0.0 ));
	
	double sourceRow, sourceCol;
	double fracRow, fracCol;
	double w1, w2, w3, w4;
	double sourceValue;
	int row, col;
	int sRow, sCol, sRowp1, sColp1;
	BOOL bSkipOutOfBounds;
	
	for ( row=0; row<m_cRows; ++row )
	{
		for ( col=0; col<m_cCols; ++col )
		{
			// the pixel at sourceRow, sourceCol is an "phantom" pixel that doesn't really exist, and
			// whose value must be manufactured from surrounding real pixels (i.e., since 
			// sourceRow and sourceCol are floating point, not ints, there's not a real pixel there)
			// The idea is that if we can calculate the value of this phantom pixel, then its 
			// displacement will exactly fit into the current pixel at row, col (which are both ints)
			
			sourceRow = (double)row - At( m_DispV, row, col );
			sourceCol = (double)col - At( m_DispH, row, col );
			
			// weights for bi-linear interpolation
			
			fracRow = sourceRow - (int)sourceRow;
			fracCol = sourceCol - (int)sourceCol;
			
			
			w1 = ( 1.0 - fracRow ) * ( 1.0 - fracCol );
			w2 = ( 1.0 - fracRow ) * fracCol;
			w3 = fracRow * ( 1 - fracCol );
			w4 = fracRow * fracCol;
			
			
			// limit indexes

/*
			while (sourceRow >= m_cRows ) sourceRow -= m_cRows;
			while (sourceRow < 0 ) sourceRow += m_cRows;
			
			while (sourceCol >= m_cCols ) sourceCol -= m_cCols;
			while (sourceCol < 0 ) sourceCol += m_cCols;
*/
			bSkipOutOfBounds = FALSE;

			if ( (sourceRow + 1.0) >= m_cRows )	bSkipOutOfBounds = TRUE;
			if ( sourceRow < 0 )				bSkipOutOfBounds = TRUE;
			
			if ( (sourceCol + 1.0) >= m_cCols )	bSkipOutOfBounds = TRUE;
			if ( sourceCol < 0 )				bSkipOutOfBounds = TRUE;
			
			if ( bSkipOutOfBounds == FALSE )
			{
				// the supporting pixels for the "phantom" source pixel are all within the 
				// bounds of the character grid.
				// Manufacture its value by bi-linear interpolation of surrounding pixels
				
				sRow = (int)sourceRow;
				sCol = (int)sourceCol;
				
				sRowp1 = sRow + 1;
				sColp1 = sCol + 1;
				
				while (sRowp1 >= m_cRows ) sRowp1 -= m_cRows;
				while (sRowp1 < 0 ) sRowp1 += m_cRows;
				
				while (sColp1 >= m_cCols ) sColp1 -= m_cCols;
				while (sColp1 < 0 ) sColp1 += m_cCols;
				
				// perform bi-linear interpolation
				
				sourceValue =	w1 * At( inputVector, sRow  , sCol   ) +
					w2 * At( inputVector, sRow  , sColp1 ) +
					w3 * At( inputVector, sRowp1, sCol   ) +
					w4 * At( inputVector, sRowp1, sColp1 );
			}
			else
			{
				// At least one supporting pixel for the "phantom" pixel is outside the
				// bounds of the character grid. Set its value to "background"

				sourceValue = 1.0;  // "background" color in the -1 -> +1 range of inputVector
			}
			
			mappedVector[ row ][ col ] = 0.5 * ( 1.0 - sourceValue );  // conversion to 0->1 range we are using for mappedVector
			
		}
	}
	
	// now, invert again while copying back into original vector
	
	for ( row=0; row<m_cRows; ++row )
	{
		for ( col=0; col<m_cCols; ++col )
		{
			At( inputVector, row, col ) = 1.0 - 2.0 * mappedVector[ row ][ col ];
		}
	}			
	
}



void CMNistDoc::CalculateNeuralNet(double *inputVector, int count, 
								   double* outputVector /* =NULL */, int oCount /* =0 */,
								   std::vector< std::vector< double > >* pNeuronOutputs /* =NULL */,
								   BOOL bDistort /* =FALSE */ )
{
	// wrapper function for neural net's Calculate() function, needed because the NN is a protected member
	// waits on the neural net mutex (using the CAutoMutex object, which automatically releases the
	// mutex when it goes out of scope) so as to restrict access to one thread at a time
	
	CAutoMutex tlo( m_utxNeuralNet );
	
	if ( bDistort != FALSE )
	{
		GenerateDistortionMap();
		ApplyDistortionMap( inputVector );
	}
	
	
	m_NN.Calculate( inputVector, count, outputVector, oCount, pNeuronOutputs );
	
}


void CMNistDoc::BackpropagateNeuralNet(double *inputVector, int iCount, double* targetOutputVector, 
									   double* actualOutputVector, int oCount, 
									   std::vector< std::vector< double > >* pMemorizedNeuronOutputs, 
									   BOOL bDistort )
{
	// function to backpropagate through the neural net. 
	
	ASSERT( (inputVector != NULL) && (targetOutputVector != NULL) && (actualOutputVector != NULL) );

	///////////////////////////////////////////////////////////////////////
	//
	// CODE REVIEW NEEDED:
	//
	// It does not seem worthwhile to backpropagate an error that's very small.  "Small" needs to be defined
	// and for now, "small" is set to a fixed size of pattern error ErrP <= 0.10 * MSE, then there will
	// not be a backpropagation of the error.  The current MSE is updated from the neural net dialog CDlgNeuralNet

	BOOL bWorthwhileToBackpropagate;  /////// part of code review
	
	{	
		// local scope for capture of the neural net, only during the forward calculation step,
		// i.e., we release neural net for other threads after the forward calculation, and after we
		// have stored the outputs of each neuron, which are needed for the backpropagation step
		
		CAutoMutex tlo( m_utxNeuralNet );
		
		// determine if it's time to adjust the learning rate
		
		if ( (( m_cBackprops % m_nAfterEveryNBackprops ) == 0) && (m_cBackprops != 0) )
		{
			double eta = m_NN.m_etaLearningRate;
			eta *= m_dEtaDecay;
			if ( eta < m_dMinimumEta )
				eta = m_dMinimumEta;
			m_NN.m_etaLearningRatePrevious = m_NN.m_etaLearningRate;
			m_NN.m_etaLearningRate = eta;
		}
		
		
		// determine if it's time to adjust the Hessian (currently once per epoch)
		
		if ( (m_bNeedHessian != FALSE) || (( m_cBackprops % ::GetPreferences().m_nItemsTrainingImages ) == 0) )
		{
			// adjust the Hessian.  This is a lengthy operation, since it must process approx 500 labels
			
			CalculateHessian();
			
			m_bNeedHessian = FALSE;
		}
		
		
		// determine if it's time to randomize the sequence of training patterns (currently once per epoch)
		
		if ( ( m_cBackprops % ::GetPreferences().m_nItemsTrainingImages ) == 0 )
		{
			RandomizeTrainingPatternSequence();
		}
		
		
		// increment counter for tracking number of backprops
		
		m_cBackprops++;
		
		
		// forward calculate through the neural net
		
		CalculateNeuralNet( inputVector, iCount, actualOutputVector, oCount, pMemorizedNeuronOutputs, bDistort );


		// calculate error in the output of the neural net
		// note that this code duplicates that found in many other places, and it's probably sensible to 
		// define a (global/static ??) function for it

		double dMSE = 0.0;
		for ( int ii=0; ii<10; ++ii )
		{
			dMSE += ( actualOutputVector[ii]-targetOutputVector[ii] ) * ( actualOutputVector[ii]-targetOutputVector[ii] );
		}
		dMSE /= 2.0;

		if ( dMSE <= ( 0.10 * m_dEstimatedCurrentMSE ) )
		{
			bWorthwhileToBackpropagate = FALSE;
		}
		else
		{
			bWorthwhileToBackpropagate = TRUE;
		}

		
		if ( (bWorthwhileToBackpropagate != FALSE) && (pMemorizedNeuronOutputs == NULL) )
		{
			// the caller has not provided a place to store neuron outputs, so we need to
			// backpropagate now, while the neural net is still captured.  Otherwise, another thread
			// might come along and call CalculateNeuralNet(), which would entirely change the neuron
			// outputs and thereby inject errors into backpropagation 
			
			m_NN.Backpropagate( actualOutputVector, targetOutputVector, oCount, NULL );
			
			SetModifiedFlag( TRUE );
			
			// we're done, so return
			
			return;
		}
		
	}
	
	// if we have reached here, then the mutex for the neural net has been released for other 
	// threads.  The caller must have provided a place to store neuron outputs, which we can 
	// use to backpropagate, even if other threads call CalculateNeuralNet() and change the outputs
	// of the neurons

	if ( (bWorthwhileToBackpropagate != FALSE) )
	{
		m_NN.Backpropagate( actualOutputVector, targetOutputVector, oCount, pMemorizedNeuronOutputs );
		
		// set modified flag to prevent closure of doc without a warning
		
		SetModifiedFlag( TRUE );
	}
	
}



void CMNistDoc::CalculateHessian()
{
	// controls the Neural network's calculation if the diagonal Hessian for the Neural net
	// This will be called from a thread, so although the calculation is lengthy, it should not interfere
	// with the UI
	
	// we need the neural net exclusively during this calculation, so grab it now
	
	CAutoMutex tlo( m_utxNeuralNet );
	
	double inputVector[841] = {0.0};  // note: 29x29, not 28x28
	double targetOutputVector[10] = {0.0};
	double actualOutputVector[10] = {0.0};
	
	unsigned char grayLevels[g_cImageSize * g_cImageSize] = { 0 };
	int label = 0;
	int ii, jj;
	UINT kk;
	
	// calculate the diagonal Hessian using 500 random patterns, per Yann LeCun 1998 "Gradient-Based Learning
	// Applied To Document Recognition"
	
	// message to dialog that we are commencing calculation of the Hessian
	
	if ( m_hWndForBackpropPosting != NULL )
	{
		// wParam == 4L -> related to Hessian, lParam == 1L -> commenced calculation
		::PostMessage( m_hWndForBackpropPosting, UWM_BACKPROPAGATION_NOTIFICATION, 4L, 1L );
	}
	
	
	// some of this code is similar to the BackpropagationThread() code
	
	m_NN.EraseHessianInformation();
	
	UINT numPatternsSampled = ::GetPreferences().m_nNumHessianPatterns ;
	
	for ( kk=0; kk<numPatternsSampled; ++kk )
	{
		GetRandomTrainingPattern( grayLevels, &label, TRUE );
		
		if ( label < 0 ) label = 0;
		if ( label > 9 ) label = 9;
		
		
		// pad to 29x29, convert to double precision
		
		for ( ii=0; ii<841; ++ii )
		{
			inputVector[ ii ] = 1.0;  // one is white, -one is black
		}
		
		// top row of inputVector is left as zero, left-most column is left as zero 
		
		for ( ii=0; ii<g_cImageSize; ++ii )
		{
			for ( jj=0; jj<g_cImageSize; ++jj )
			{
				inputVector[ 1 + jj + 29*(ii+1) ] = (double)((int)(unsigned char)grayLevels[ jj + g_cImageSize*ii ])/128.0 - 1.0;  // one is white, -one is black
			}
		}
		
		// desired output vector
		
		for ( ii=0; ii<10; ++ii )
		{
			targetOutputVector[ ii ] = -1.0;
		}
		targetOutputVector[ label ] = 1.0;
		
		
		// apply distortion map to inputVector.  It's not certain that this is needed or helpful.
		// The second derivatives do NOT rely on the output of the neural net (i.e., because the 
		// second derivative of the MSE function is exactly 1 (one), regardless of the actual output
		// of the net).  However, since the backpropagated second derivatives rely on the outputs of
		// each neuron, distortion of the pattern might reveal previously-unseen information about the
		// nature of the Hessian.  But I am reluctant to give the full distortion, so I set the
		// severityFactor to only 2/3 approx
		
		GenerateDistortionMap( 0.65 );
		ApplyDistortionMap( inputVector );
		
		
		// forward calculate the neural network
		
		m_NN.Calculate( inputVector, 841, actualOutputVector, 10, NULL );
		
		
		// backpropagate the second derivatives
		
		m_NN.BackpropagateSecondDervatives( actualOutputVector, targetOutputVector, 10 );
		
		
		// progress message to dialog that we are calculating the Hessian
		
		if ( kk%50 == 0 )
		{
			// every 50 iterations ...
			if ( m_hWndForBackpropPosting != NULL )
			{
				// wParam == 4L -> related to Hessian, lParam == 2L -> progress indicator
				::PostMessage( m_hWndForBackpropPosting, UWM_BACKPROPAGATION_NOTIFICATION, 4L, 2L );
			}
		}
		
		if ( m_bBackpropThreadAbortFlag != FALSE )
			break;
		
	}
	
	m_NN.DivideHessianInformationBy( (double)numPatternsSampled );
	
	// message to dialog that we are finished calculating the Hessian
	
	if ( m_hWndForBackpropPosting != NULL )
	{
		// wParam == 4L -> related to Hessian, lParam == 4L -> finished calculation
		::PostMessage( m_hWndForBackpropPosting, UWM_BACKPROPAGATION_NOTIFICATION, 4L, 4L );
	}
	
}



BOOL CMNistDoc::CanCloseFrame(CFrameWnd* pFrame) 
{
	// check if any threads are running before we allow the main frame to close down
	
	BOOL bRet = TRUE;
	
	if ( (m_bBackpropThreadsAreRunning != FALSE) || (m_bTestingThreadsAreRunning != FALSE) )
	{
		CString str;
		int iRet;
		
		str.Format( _T( "This will stop backpropagation and/or testing threads \n" )
			_T( "Are you certain that you wish to close the application? \n\n" )
			_T( "Click \"Yes\" to stop all threads and close the application \n" )
			_T( "Click \"No\" or \"Cancel\" to continue running the threads and the application " ) );
		
		iRet = ::MessageBox( NULL, str, _T( "Threads Are Running" ), MB_ICONEXCLAMATION | MB_YESNOCANCEL );
		
		if ( (iRet == IDYES) || (iRet == IDOK) )
		{
			bRet = TRUE;
			
			if ( m_bBackpropThreadsAreRunning != FALSE )
			{
				StopBackpropagation();
			}
			
			if ( m_bTestingThreadsAreRunning != FALSE )
			{
				StopTesting();
			}
		}
		else
		{
			bRet = FALSE;
		}
	}
	
	if ( bRet != FALSE )
	{
		// only call the base class if, so far, it's still safe to close.  If we
		// always called the base class, then we would get needless reminders to save the
		// current document.  These reminders are not needed, since we're not closing 
		// anyway
		
		bRet &= COleDocument::CanCloseFrame(pFrame);
	}
	
	return bRet;
}



BOOL CMNistDoc::StartBackpropagation(UINT iStartPattern /* =0 */, UINT iNumThreads /* =2 */, HWND hWnd /* =NULL */,
									 double initialEta /* =0.005 */, double minimumEta /* =0.000001 */, double etaDecay /* =0.990 */, 
									 UINT nAfterEvery  /* =1000 */, BOOL bDistortPatterns /* =TRUE */, double estimatedCurrentMSE /* =1.0 */)
{
	if ( m_bBackpropThreadsAreRunning != FALSE )
		return FALSE;
	
	m_bBackpropThreadAbortFlag = FALSE;
	m_bBackpropThreadsAreRunning = TRUE;
	m_iNumBackpropThreadsRunning = 0;
	m_iBackpropThreadIdentifier = 0;
	m_cBackprops = iStartPattern;
	m_bNeedHessian = TRUE;
	
	m_iNextTrainingPattern = iStartPattern;
	m_hWndForBackpropPosting = hWnd;
	
	if ( m_iNextTrainingPattern < 0 ) 
		m_iNextTrainingPattern = 0;
	if ( m_iNextTrainingPattern >= ::GetPreferences().m_nItemsTrainingImages )
		m_iNextTrainingPattern = ::GetPreferences().m_nItemsTrainingImages - 1;
	
	if ( iNumThreads < 1 ) 
		iNumThreads = 1;
	if ( iNumThreads > 10 )  // 10 is arbitrary upper limit
		iNumThreads = 10;
	
	m_NN.m_etaLearningRate = initialEta;
	m_NN.m_etaLearningRatePrevious = initialEta;
	m_dMinimumEta = minimumEta;
	m_dEtaDecay = etaDecay;
	m_nAfterEveryNBackprops = nAfterEvery;
	m_bDistortTrainingPatterns = bDistortPatterns;

	m_dEstimatedCurrentMSE = estimatedCurrentMSE;  // estimated number that will define whether a forward calculation's error is significant enough to warrant backpropagation

	RandomizeTrainingPatternSequence();
	
	for ( UINT ii=0; ii<iNumThreads; ++ii )
	{
		CWinThread* pThread = ::AfxBeginThread( BackpropagationThread, (LPVOID)this, 
			THREAD_PRIORITY_BELOW_NORMAL, 0, CREATE_SUSPENDED, NULL );
		
		if ( pThread == NULL )
		{
			// creation failed; un-do everything
			
			StopBackpropagation();
			return FALSE;
		}
		
		pThread->m_bAutoDelete = FALSE;		 
		m_pBackpropThreads[ ii ] = pThread;
		pThread->ResumeThread();
		m_iNumBackpropThreadsRunning++;
	}
	
	return TRUE;
	
}




void CMNistDoc::StopBackpropagation()
{
	// stops all the backpropagation threads
	
	if ( m_bBackpropThreadsAreRunning == FALSE )
	{
		// it's curious to select "stop" if no threads are running, but perform some
		// shutdown safeguards, just to be certain
		
		m_bBackpropThreadAbortFlag = TRUE;
		m_bBackpropThreadsAreRunning = FALSE;
		m_iNumBackpropThreadsRunning = 0;
		m_iBackpropThreadIdentifier = 0;
		m_cBackprops = 0;
		
		return;
	}
	
	m_bBackpropThreadAbortFlag = TRUE;
	
	HANDLE hThread[25];
	CString msg;
	DWORD dwTimeOut = 5000;  // 5 second default timeout
	UINT ii;
	
	while ( m_iNumBackpropThreadsRunning > 0 )
	{
		for ( ii=0; ii<m_iNumBackpropThreadsRunning; ++ii )
		{
			hThread[ ii ] = m_pBackpropThreads[ ii ]->m_hThread;
		}
		
		DWORD dwRet = ::WaitForMultipleObjects( m_iNumBackpropThreadsRunning, hThread, FALSE, dwTimeOut );
		DWORD dwErr = ::GetLastError();		// if an error occurred get its value as soon as possible
		
		if ( dwRet==WAIT_FAILED )
		{
			LPVOID lpMsgBuf;
			::FormatMessage( FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM | FORMAT_MESSAGE_IGNORE_INSERTS,
				NULL, dwErr, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT), (LPTSTR) &lpMsgBuf, 0, NULL );
			
			::MessageBox( NULL, (LPCTSTR)lpMsgBuf, _T("Error Waiting For Backpropagation Thread Shutdown"), MB_OK | MB_ICONINFORMATION );
			LocalFree( lpMsgBuf );
		}
		else if ( dwRet==WAIT_TIMEOUT )
		{
			// bad -- no threads are responding
			// give user option of waiting a bit more, or of terminating the threads forcefully
			
			msg.Format( _T("No thread has responded after waiting %d milliseconds\n\n")
				_T("Click \"Retry\" to wait another %d milliseconds and give a thread\n")
				_T("a chance to respond\n\n")
				_T("Click \"Cancel\" to terminate a thread forcibly\n\n")
				_T("Note: Forceful termination is not recommended and might cause memory leaks"), dwTimeOut, dwTimeOut );
			if ( IDCANCEL == ::MessageBox( NULL, msg, _T("No Thread Is Responding"), MB_RETRYCANCEL|MB_ICONEXCLAMATION|MB_DEFBUTTON1 ) )
			{
				// forceful thread termination was selected
				// pick first thread from list and terminate it
				
				::TerminateThread( hThread[0], 98 );	// specify exit code of "98"
			}
			
		}
		else if ( dwRet>=WAIT_ABANDONED_0 && dwRet<=(WAIT_ABANDONED_0+m_iNumBackpropThreadsRunning-1) )
		{
			msg.Format( _T("Thread reports mutex was abandoned") );
			::MessageBox( NULL, msg, _T("Mutex Abandoned"), MB_OK|MB_ICONEXCLAMATION );
		}
		else if ( dwRet>=WAIT_OBJECT_0 && dwRet<=(WAIT_OBJECT_0+m_iNumBackpropThreadsRunning-1) )
		{
			// the most common and expected return value
			// delete the object that signalled
			
			int nDex = dwRet - WAIT_OBJECT_0;
			
			delete m_pBackpropThreads[ nDex ];
			
			for ( ii=nDex; ii<m_iNumBackpropThreadsRunning-1; ++ii )
			{
				m_pBackpropThreads[ ii ] = m_pBackpropThreads[ ii+1 ];
			}
			
			m_iNumBackpropThreadsRunning--;
			
		}
		else
		{
			ASSERT( FALSE );	// shouldn't be able to get here
		}
	}
	
	
	// at this point, all threads have been terminated, so re-set flags to allow
	// for future re-start of the threads
	
	m_bBackpropThreadAbortFlag = TRUE;
	m_bBackpropThreadsAreRunning = FALSE;
	m_iNumBackpropThreadsRunning = 0;
	m_iBackpropThreadIdentifier = 0;
	m_cBackprops = 0;
	
	
	
}




UINT CMNistDoc::BackpropagationThread(LPVOID pVoid)
{
	// thread for backpropagation training of NN
	//
	// thread is "owned" by the doc, and accepts a pointer to the doc
	// continuously backpropagates until m_bThreadAbortFlag is set to TRUE
	
	CMNistDoc* pThis = reinterpret_cast< CMNistDoc* >( pVoid );
	
	ASSERT( pThis != NULL );
	
	// set thread name (helps during debugging)
	
	char str[25] = {0};  // must use chars, not TCHARs, for SetThreadname function
	sprintf( str, "BACKP%02d", pThis->m_iBackpropThreadIdentifier++ );
	SetThreadName( -1, str );
	
	// do the work
	
	double inputVector[841] = {0.0};  // note: 29x29, not 28x28
	double targetOutputVector[10] = {0.0};
	double actualOutputVector[10] = {0.0};
	double dMSE;
	UINT scaledMSE;
	unsigned char grayLevels[g_cImageSize * g_cImageSize] = { 0 };
	int label = 0;
	int ii, jj;
	UINT iSequentialNum;
	
	std::vector< std::vector< double > > memorizedNeuronOutputs;
	
	while ( pThis->m_bBackpropThreadAbortFlag == FALSE )
	{
		int iRet = pThis->GetNextTrainingPattern( grayLevels, &label, TRUE, TRUE, &iSequentialNum );
		
		if ( label < 0 ) label = 0;
		if ( label > 9 ) label = 9;
		
		// post message to the dialog, telling it which pattern this thread is currently working on
		
		if ( pThis->m_hWndForBackpropPosting != NULL )
		{
			::PostMessage( pThis->m_hWndForBackpropPosting, UWM_BACKPROPAGATION_NOTIFICATION, 1L, (LPARAM)iSequentialNum );
		}
		
		
		// pad to 29x29, convert to double precision
		
		for ( ii=0; ii<841; ++ii )
		{
			inputVector[ ii ] = 1.0;  // one is white, -one is black
		}
		
		// top row of inputVector is left as zero, left-most column is left as zero 
		
		for ( ii=0; ii<g_cImageSize; ++ii )
		{
			for ( jj=0; jj<g_cImageSize; ++jj )
			{
				inputVector[ 1 + jj + 29*(ii+1) ] = (double)((int)(unsigned char)grayLevels[ jj + g_cImageSize*ii ])/128.0 - 1.0;  // one is white, -one is black
			}
		}
		
		// desired output vector
		
		for ( ii=0; ii<10; ++ii )
		{
			targetOutputVector[ ii ] = -1.0;
		}
		targetOutputVector[ label ] = 1.0;
		
		// now backpropagate

		pThis->BackpropagateNeuralNet( inputVector, 841, targetOutputVector, actualOutputVector, 10, 
			&memorizedNeuronOutputs, pThis->m_bDistortTrainingPatterns );

		
		// calculate error for this pattern and post it to the hwnd so it can calculate a running 
		// estimate of MSE
		
		dMSE = 0.0;
		for ( ii=0; ii<10; ++ii )
		{
			dMSE += ( actualOutputVector[ii]-targetOutputVector[ii] ) * ( actualOutputVector[ii]-targetOutputVector[ii] );
		}
		dMSE /= 2.0;

		scaledMSE = (UINT)( sqrt( dMSE ) * 2.0e8 );  // arbitrary large pre-agreed upon scale factor; taking sqrt is simply to improve the scaling
		
		if ( pThis->m_hWndForBackpropPosting != NULL )
		{
			::PostMessage( pThis->m_hWndForBackpropPosting, UWM_BACKPROPAGATION_NOTIFICATION, 2L, (LPARAM)scaledMSE );
		}
		
		
		// determine the neural network's answer, and compare it to the actual answer.
		// Post a message if the answer was incorrect, so the dialog can display mis-recognition
		// statistics
		
		int iBestIndex = 0;
		double maxValue = -99.0;
		
		for ( ii=0; ii<10; ++ii )
		{
			if ( actualOutputVector[ ii ] > maxValue )
			{
				iBestIndex = ii;
				maxValue = actualOutputVector[ ii ];
			}
		}

		if ( iBestIndex != label )
		{
			// pattern was mis-recognized.  Notify the testing dialog
			if ( pThis->m_hWndForBackpropPosting != NULL )
			{
				::PostMessage( pThis->m_hWndForBackpropPosting, UWM_BACKPROPAGATION_NOTIFICATION, 8L, (LPARAM)0L );
			}
		}
		
	}  // end of main "while not abort flag" loop
	
	return 0L;
}







BOOL CMNistDoc::StartTesting(UINT iStartingPattern, UINT iNumThreads, HWND hWnd, BOOL bDistortPatterns,
							 UINT iWhichImageSet /* =1 */ )
{
	// creates and starts testing threads
	
	if ( m_bTestingThreadsAreRunning != FALSE )
		return FALSE;
	
	m_bTestingThreadAbortFlag = FALSE;
	m_bTestingThreadsAreRunning = TRUE;
	m_iNumTestingThreadsRunning = 0;
	m_iTestingThreadIdentifier = 0;
	
	m_iNextTestingPattern = iStartingPattern;
	m_hWndForTestingPosting = hWnd;
	m_iWhichImageSet = iWhichImageSet;
	
	if ( m_iWhichImageSet > 1 )
		m_iWhichImageSet = 1;
	if ( m_iWhichImageSet < 0 )  // which is not possible, since m_iWhichImageSet is a UINT
		m_iWhichImageSet = 0;
	
	if ( m_iNextTestingPattern < 0 ) 
		m_iNextTestingPattern = 0;
	if ( m_iNextTestingPattern >= ::GetPreferences().m_nItemsTestingImages )
		m_iNextTestingPattern = ::GetPreferences().m_nItemsTestingImages - 1;
	
	if ( iNumThreads < 1 ) 
		iNumThreads = 1;
	if ( iNumThreads > 10 )  // 10 is arbitrary upper limit
		iNumThreads = 10;
	
	m_bDistortTestingPatterns = bDistortPatterns;
	
	for ( UINT ii=0; ii<iNumThreads; ++ii )
	{
		CWinThread* pThread = ::AfxBeginThread( TestingThread, (LPVOID)this, 
			THREAD_PRIORITY_BELOW_NORMAL, 0, CREATE_SUSPENDED, NULL );
		
		if ( pThread == NULL )
		{
			// creation failed; un-do everything
			
			StopTesting();
			return FALSE;
		}
		
		pThread->m_bAutoDelete = FALSE;		 
		m_pTestingThreads[ ii ] = pThread;
		pThread->ResumeThread();
		m_iNumTestingThreadsRunning++;
	}
	
	return TRUE;
}

void CMNistDoc::StopTesting()
{
	// stops all the testing threads
	
	if ( m_bTestingThreadsAreRunning == FALSE )
	{
		// it's curious to select "stop" if no threads are running, but perform some
		// shutdown safeguards, just to be certain
		
		m_bTestingThreadAbortFlag = TRUE;
		m_bTestingThreadsAreRunning = FALSE;
		m_iNumTestingThreadsRunning = 0;
		m_iTestingThreadIdentifier = 0;
		return;
	}
	
	m_bTestingThreadAbortFlag = TRUE;
	
	HANDLE hThread[25];
	CString msg;
	DWORD dwTimeOut = 5000;  // 5 second default timeout
	UINT ii;
	
	while ( m_iNumTestingThreadsRunning > 0 )
	{
		for ( ii=0; ii<m_iNumTestingThreadsRunning; ++ii )
		{
			hThread[ ii ] = m_pTestingThreads[ ii ]->m_hThread;
		}
		
		DWORD dwRet = ::WaitForMultipleObjects( m_iNumTestingThreadsRunning, hThread, FALSE, dwTimeOut );
		DWORD dwErr = ::GetLastError();		// if an error occurred get its value as soon as possible
		
		if ( dwRet==WAIT_FAILED )
		{
			LPVOID lpMsgBuf;
			::FormatMessage( FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM | FORMAT_MESSAGE_IGNORE_INSERTS,
				NULL, dwErr, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT), (LPTSTR) &lpMsgBuf, 0, NULL );
			
			::MessageBox( NULL, (LPCTSTR)lpMsgBuf, _T("Error Waiting For Testing Thread Shutdown"), MB_OK | MB_ICONINFORMATION );
			LocalFree( lpMsgBuf );
		}
		else if ( dwRet==WAIT_TIMEOUT )
		{
			// bad -- no threads are responding
			// give user option of waiting a bit more, or of terminating the threads forcefully
			
			msg.Format( _T("No thread has responded after waiting %d milliseconds\n\n")
				_T("Click \"Retry\" to wait another %d milliseconds and give a thread\n")
				_T("a chance to respond\n\n")
				_T("Click \"Cancel\" to terminate a thread forcibly\n\n")
				_T("Note: Forceful termination is not recommended and might cause memory leaks"), dwTimeOut, dwTimeOut );
			if ( IDCANCEL == ::MessageBox( NULL, msg, _T("No Thread Is Responding"), MB_RETRYCANCEL|MB_ICONEXCLAMATION|MB_DEFBUTTON1 ) )
			{
				// forceful thread termination was selected
				// pick first thread from list and terminate it
				
				::TerminateThread( hThread[0], 98 );	// specify exit code of "98"
			}
			
		}
		else if ( dwRet>=WAIT_ABANDONED_0 && dwRet<=(WAIT_ABANDONED_0+m_iNumTestingThreadsRunning-1) )
		{
			msg.Format( _T("Thread reports mutex was abandoned") );
			::MessageBox( NULL, msg, _T("Mutex Abandoned"), MB_OK|MB_ICONEXCLAMATION );
		}
		else if ( dwRet>=WAIT_OBJECT_0 && dwRet<=(WAIT_OBJECT_0+m_iNumTestingThreadsRunning-1) )
		{
			// the most common and expected return value
			// delete the object that signalled
			
			int nDex = dwRet - WAIT_OBJECT_0;
			
			delete m_pTestingThreads[ nDex ];
			
			for ( ii=nDex; ii<m_iNumTestingThreadsRunning-1; ++ii )
			{
				m_pTestingThreads[ ii ] = m_pTestingThreads[ ii+1 ];
			}
			
			m_iNumTestingThreadsRunning--;
			
		}
		else
		{
			ASSERT( FALSE );	// shouldn't be able to get here
		}
	}
	
	
	// at this point, all threads have been terminated, so re-set flags to allow
	// for future re-start of the threads
	
	m_bTestingThreadAbortFlag = TRUE;
	m_bTestingThreadsAreRunning = FALSE;
	m_iNumTestingThreadsRunning = 0;
	m_iTestingThreadIdentifier = 0;
	
	
}

UINT CMNistDoc::TestingThread(LPVOID pVoid)
{
	// thread for testing of Neural net
	// Continuously get the doc's next pattern, puts it through the neural net, and
	// inspects the output.  As the thread goes through the patterns, it post messages to the
	// m_hWndForTestingPosting, which presumably is the dialog that shows testing results,
	// advising it of the current pattern being tested.  If the actual output from the 
	// neural net differs from the desired output, another message is posted, advising the 
	// m_hWndForTestingPosting of the identity of the mis-recognized pattern
	
	// thread is owned by the doc and accepts a pointer to the doc as a parameter
	
	
	CMNistDoc* pThis = reinterpret_cast< CMNistDoc* >( pVoid );
	
	ASSERT( pThis != NULL );
	
	// set thread name (helps during debugging)
	
	char str[25] = {0};  // must use chars, not TCHARs, for SetThreadname function
	sprintf( str, "TEST%02d", pThis->m_iTestingThreadIdentifier++ );
	SetThreadName( -1, str );
	
	// do the work
	
	double inputVector[841] = {0.0};  // note: 29x29, not 28x28
	double targetOutputVector[10] = {0.0};
	double actualOutputVector[10] = {0.0};

	double dPatternMSE = 0.0;
	double dTotalMSE = 0.0;
	UINT scaledMSE = 0;
	UINT iPatternsProcessed = 0;
	
	unsigned char grayLevels[g_cImageSize * g_cImageSize] = { 0 };
	int label = 0;
	int ii, jj;
	UINT iPatNum, iSequentialNum;
	
	while ( pThis->m_bTestingThreadAbortFlag == FALSE )
	{
		// testing image set or training image set
		
		if ( pThis->m_iWhichImageSet == 1 )
		{
			// testing set
			
			iPatNum = pThis->GetNextTestingPattern( grayLevels, &label, TRUE );
			
			// post message to the dialog, telling it which pattern this thread is currently working on
			
			if ( pThis->m_hWndForTestingPosting != NULL )
			{
				::PostMessage( pThis->m_hWndForTestingPosting, UWM_TESTING_NOTIFICATION, 1L, (LPARAM)iPatNum );
			}
		}
		else
		{
			// training set
			
			iPatNum = pThis->GetNextTrainingPattern( grayLevels, &label, TRUE, FALSE, &iSequentialNum );
			
			// post message to the dialog, telling it which pattern this thread is currently working on
			
			if ( pThis->m_hWndForTestingPosting != NULL )
			{
				::PostMessage( pThis->m_hWndForTestingPosting, UWM_TESTING_NOTIFICATION, 1L, (LPARAM)iSequentialNum );
			}
		}
		
		
		if ( label < 0 ) label = 0;
		if ( label > 9 ) label = 9;
		
		
		
		
		// pad to 29x29, convert to double precision
		
		for ( ii=0; ii<841; ++ii )
		{
			inputVector[ ii ] = 1.0;  // one is white, -one is black
		}
		
		// top row of inputVector is left as zero, left-most column is left as zero 
		
		for ( ii=0; ii<g_cImageSize; ++ii )
		{
			for ( jj=0; jj<g_cImageSize; ++jj )
			{
				inputVector[ 1 + jj + 29*(ii+1) ] = (double)((int)(unsigned char)grayLevels[ jj + g_cImageSize*ii ])/128.0 - 1.0;  // one is white, -one is black
			}
		}
		
		// desired output vector
		
		for ( ii=0; ii<10; ++ii )
		{
			targetOutputVector[ ii ] = -1.0;
		}
		targetOutputVector[ label ] = 1.0;
		

		// now calculate output of neural network
		
		pThis->CalculateNeuralNet( inputVector, 841, actualOutputVector, 10, NULL, pThis->m_bDistortTestingPatterns );

		
		// calculate error for this pattern and accumulate it for posting of
		// total MSE of all patterns when thread is exiting
		
		dPatternMSE = 0.0;
		for ( ii=0; ii<10; ++ii )
		{
			dPatternMSE += ( actualOutputVector[ii]-targetOutputVector[ii] ) * ( actualOutputVector[ii]-targetOutputVector[ii] );
		}
		dPatternMSE /= 2.0;

		dTotalMSE += dPatternMSE;
		++iPatternsProcessed;		
		

		// determine the neural network's answer, and compare it to the actual answer
		
		int iBestIndex = 0;
		double maxValue = -99.0;
		UINT code;
		
		for ( ii=0; ii<10; ++ii )
		{
			if ( actualOutputVector[ ii ] > maxValue )
			{
				iBestIndex = ii;
				maxValue = actualOutputVector[ ii ];
			}
		}
		
		
		// moment of truth: Did neural net get the correct answer
		
		if ( iBestIndex != label )
		{
			// pattern was mis-recognized.  Notify the testing dialog
			
			// lParam is built to contain a coded bit pattern, as follows:
			//
			//  0          1          2         3
			//  0123456 7890123 456789012345678901
			// |  act  |  tar  |    pattern num   |
			//
			// where act == actual output of the neural net, and tar == target
			// this gives 2^7 = 128 possible outputs (only 10 are needed here... future expansion??)
			// and 2^18 = 262144 possible pattern numbers ( only 10000 are needed here )
			
			code  = ( iPatNum    & 0x0003FFFF );
			code |= ( label      & 0x0000007F ) << 18;
			code |= ( iBestIndex & 0x0000007F ) << 25;
			
			if ( pThis->m_hWndForTestingPosting != NULL )
			{
				::PostMessage( pThis->m_hWndForTestingPosting, UWM_TESTING_NOTIFICATION, 2L, (LPARAM)code );
			}
		}
		
	}


	// post the total MSE of tested patterns to the hwnd

	double divisor = (double)( (iPatternsProcessed>1) ? iPatternsProcessed : 1 );
	dTotalMSE /= divisor;
	scaledMSE = (UINT)( sqrt( dTotalMSE ) * 2.0e8 );  // arbitrary large pre-agreed upon scale factor; taking sqrt is simply to improve the scaling
		
	if ( pThis->m_hWndForTestingPosting != NULL )
	{
		::PostMessage( pThis->m_hWndForTestingPosting, UWM_TESTING_NOTIFICATION, 4L, (LPARAM)scaledMSE );
	}
	
	return 0L;
	
}


#undef UNIFORM_ZERO_THRU_ONE

By viewing downloads associated with this article you agree to the Terms of Service and the article's licence.

If a file you wish to view isn't highlighted, and is a text file (not binary), please let us know and we'll add colourisation support for it.

License

This article has no explicit license attached to it but may contain usage terms in the article text or the download files themselves. If in doubt please contact the author via the discussion board below.

A list of licenses authors might use can be found here

Share

About the Author

Mike O'Neill

United States United States
Mike O'Neill is a patent attorney in Southern California, where he specializes in computer and software-related patents. He programs as a hobby, and in a vain attempt to keep up with and understand the technology of his clients.

| Advertise | Privacy | Mobile
Web01 | 2.8.140827.1 | Last Updated 5 Dec 2006
Article Copyright 2006 by Mike O'Neill
Everything else Copyright © CodeProject, 1999-2014
Terms of Service
Layout: fixed | fluid