/************************************
 * Simulation of Pike's Memory Model *
 *************************************
 *
 * October 1992
 *
 * David Diller,
 * David Huber & Peter Nobel
 * Indiana University
 * Bloomington, IN  47405
 *
 * ddiller@silver.indiana.edu,
 * dhuber@ucs.indiana.edu & pnobel@ucs.indiana.edu
 *
 * This program simulates the following test-conditions:
 *      single item recognition
 *      associative recognition: intact - rearranged
 *      associative recognition: intact - mixed
 *      associative recognition: intact - new
 *      cued recall
 *
 */
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <sys/types.h>   /* these headers aren't needed when using turbo c */
#include <sys/timeb.h>   /* as above */
#include <malloc.h>                              /* for turbo c use alloc.h */
#include <float.h>

#define MAX_ENTRIES 10000

typedef struct {
    char *name;
    char *val;
} entry;

char *makeword(char *line, char stop);
char *fmakeword(FILE *f, char stop, int *len);
char x2c(char *what);
void unescape_url(char *url);
void plustospace(char *str);



#define RAND_MAX 2147483648.0         /* the maximum number returned by rand */
#define NSIM 500                           /* maximum number of simulations */
#define NMAX 32                        /* maximum number of pairs on a list */
#define SMAXS 5                            /* maximum strength strong items */
#define SMAXW 3                              /* maximum strength weak items */

#define TRUE 1

#define integer int                             /* used in d-prime function */
#define real float
static integer c__9 = 9;
static integer c__1 = 1;
static integer c__4 = 4;
float hrate;
float farate;
float dp;

FILE *fp_out,
		 *fp;                                    /* file pointer to output file */
char fout[100];                                         /* output file name */

/***********************
 * function prototypes *
 ***********************/
void menu(void);
void making_matrix(void);
void recognition_single(void);
void recognition_ass_rearr(void);
void recognition_ass_mixed(void);
void recognition_ass_new(void);
void cued_recall(void);
void ranset(void);
void perm(int *, int);
void simulate(void);
void results(void);
void display(char **, int, int);
void settings(void);
void screen(void);
void move(int, int);
void action1(int);
void action2(int);
void action3(int);
void action4(int);
void action5(int);
void action6(int);
void study(void);
float rand1(void);
int randw(int);
int getcode(void);
float single_fam1(int);
float single_fam2(int);
float double_fam(int, int);
float dot_product ();
double r_sign();
float gau(double, double);
void calculate_criterion ();

/********************
 * menu definitions *
 ********************/
char *menu1[] =
			{ "Select type of test",
				"Set length of list",
				"Set strength of list",
				"Set parameters",
				"Simulation",
				"Write results",
				"Exit",                   };

char *menu2[] =
			{ "Single item recognition",
				"Assoc. recog. - rearranged",
				"Assoc. recog. - mixed",
				"Assoc. recog. - new",
				"Cued recall",
				"Return to Main Menu",    };

char *menu3[] =
			{ "Set length of list (number of pairs)",
				"Return to Main Menu",    };

char *menu4[] =
			{ "Pure weak",
				"Pure strong",
				"Mixed",
				"Return to Main Menu",    };

char *menu5[] =
			{ "Criterion (recog.)",
        "Criterion (recall)",
				"Number of Features",
				"Storage Probability",
				"Element Mean",
				"Element Power",
				"Return to Main Menu",    };

char *menu6[] =
			{ "Number of simulations",
							"Simulate",
							"Return to Main Menu",    };

/***************************************
 * global variables for the simulation *
 ***************************************/
int nsim=50;                                       /* number of simulations */
int cond=4;                                                 /* type of test */
int n=8;             /* the length of the list (# of pairs): multiples of 8 */
int length;                               /* actual number of presentations */
int strength=1; /*strength of the list: 1-pure weak, 2-pure strong, 3-mixed */
int s1=1;                                         /* strength of weak items */
int s2=2;                                       /* strength of strong items */

/**********************************
 * global variables for the model *
 **********************************/

float items[64][64];      /* the lexicon could consist of as many as 64 two
										 item pairs with each item containing 32 features.
										 these items are randomly determined at the start and
										 stay static for checking purposes throughout. these
										 items are 'clean' and so the features are -1 or 1 */


float new_memory[32][32];  /* Matrix containing the traces stored in memory
													 There is a maxmimum of 32 possible features for
													 every item in memory. */

float matrix1[32][32];
float matrix2[32][32];


int count,                         /* counter for the number of simulations */
		count_cw,                   /* counter for weak corrects in cued recall */
		count_cs,                 /* counter for strong corrects in cued recall */
		count_hw,                       /* counter for weak hits in recognition */
		count_hs,                     /* counter for strong hits in recognition */
		count_f,                     /* counter for false alarms in recognition */
		count_fs;                      /* counter for strong false alarms (ass) */

float famtw,                         /* total sum familiarity: weak targets */
			famtw2,                /* total sum familiarity squared: weak targets */
			famts,                       /* total sum familiarity: strong targets */
			famts2,              /* total sum familiarity squared: strong targets */
			famd,                           /* total sum familiarity: distractors */
			famd2,                  /* total sum familiarity squared: distractors */
			famds,                   /* total sum familiarity: strong distractors */
			famds2;          /* total sum familiarity squared: strong distractors */

float VARTW[NSIM],          /* array for variance familiarity: weak targets */
			VARTS[NSIM],        /* array for variance familiarity: strong targets */
			VARD[NSIM],            /* array for variance familiarity: distractors */
			VARDS[NSIM],           /* array for variance familiarity: distractors */
			DPRPW[NSIM],         /* array for d-prime (weak) based on proportions */
			DPRPS[NSIM],       /* array for d-prime (strong) based on proportions */
			DPRDW[NSIM],       /* array for d-prime (weak) based on distributions */
			DPRDS[NSIM];     /* array for d-prime (strong) based on distributions */


/****************************
 * parameters for the model *
 ****************************/

#define ON 0
#define OFF 1
#define NONE -1.0

int   features=16;        /* number of features chosen */
float learn=.9;           /* probability of storing a feature */
float crit=0.45;          /* cut-off criterion used for recognition */
float recall_crit=0.0;    /* cut-off criterion used for recall */
float corw;               /* output variable of average correlation for weak targets */
float cors;               /* output variable of average correlation for strong targets */
float mycorw, mycors;
float mu     = 0.0,
	  sigma  = 1.0;

/************************
 * variables for output *
 ************************/

/* cued recall */
float pcw;                               /* final result recall: weak items */
float pcs;                             /* final result recall: strong items */

/* recognition: variance and d' for both averaged and total */
float hw;                            /* final result recognition: weak hits */
float mtw;                                 /* mean familiarity weak targets */
float vartw;                    /* variance familiarity weak targets (avg.) */
float tvartw;                  /* variance familiarity weak targets (total) */
float hs;                          /* final result recognition: strong hits */
float mts;                               /* mean familiarity strong targets */
float varts;                  /* variance familiarity strong targets (avg.) */
float tvarts;                /* variance familiarity strong targets (total) */
float f;                          /* final result recognition: false alarms */
float fs;                  /* final result recognition: strong false alarms */
float md;                                   /* mean familiarity distractors */
float mds;                           /* mean familiarity strong distractors */
float vard;                      /* variance familiarity distractors (avg.) */
float vards;              /* variance familiarity strong distractors (avg.) */
float tvard;                    /* variance familiarity distractors (total) */
float tvards;            /* variance familiarity strong distractors (total) */
float dprpw;              /* dprime based on proportions: weak items (avg.) */
float tdprpw;            /* dprime based on proportions: weak items (total) */
float dprps;            /* dprime based on proportions: strong items (avg.) */
float tdprps;          /* dprime based on proportions: strong items (total) */
float dprdw;            /* dprime based on distributions: weak items (avg.) */
float tdprdw;          /* dprime based on distributions: weak items (total) */
float dprds;          /* dprime based on distributions: strong items (avg.) */
float tdprds;        /* dprime based on distributions: strong items (total) */

/********
 * MAIN *
 ********/
main(int argc, char *argv[]) {

    entry entries[MAX_ENTRIES];
    register int x,m=0;
    int cl;

    printf("Content-type: text/html%c%c",10,10);

    if(strcmp(getenv("REQUEST_METHOD"),"POST")) {
        printf("This script should be referenced with a METHOD of POST.\n");
        printf("If you don't understand this, see this ");
        printf("<A HREF=\"http://www.ncsa.uiuc.edu/SDG/Software/Mosaic/Docs/fill-out-forms/overview.htm
l\">forms overview</A>.%c",10);
        exit(1);
    }
    if(strcmp(getenv("CONTENT_TYPE"),"application/x-www-form-urlencoded")) {
        printf("This script can only be used to decode form results. \n");
        exit(1);
    }
    cl = atoi(getenv("CONTENT_LENGTH"));

    for(x=0;cl && (!feof(stdin));x++) {
        m=x;
        entries[x].val = fmakeword(stdin,'&',&cl);
        plustospace(entries[x].val);
        unescape_url(entries[x].val);
        entries[x].name = makeword(entries[x].val,'=');

    if (strcmp(entries[x].name, "TypeOfTest") == 0){
        if (strcmp(entries[x].val, "Single item recognition") == 0)
            cond = 1;

        else if (strcmp(entries[x].val, "Associative recognition - rearranged") == 0)
            cond = 2;
        else if (strcmp(entries[x].val, "Associative recognition - mixed") == 0)
            cond = 3;
        else if (strcmp(entries[x].val, "Associative recognition - new") == 0)
            cond = 4;
        else if (strcmp(entries[x].val, "Cued recall") == 0)
            cond = 5;
        }
    else if (strcmp(entries[x].name, "ListLength") == 0){
        sscanf(entries[x].val, "%d", &n);
        if (n > NMAX) {
            printf("Maximum number of pairs is %d!", NMAX);
            n = NMAX;
            }
        if (n % 8 != 0) {                           /* if not a multiple of 8 */
            printf("Number of pairs should be a multiple of 8! Setting list length to 8.");
            n=8;
            }
        }
    else if (strcmp(entries[x].name, "ListStrength") == 0){
        if (strcmp(entries[x].val, "Pure Weak") == 0)
            strength = 1;
        else if (strcmp(entries[x].val, "Pure Strong") == 0)
            strength = 2;
        else if (strcmp(entries[x].val, "Mixed") == 0)
            strength = 3;
        }
    else if (strcmp(entries[x].name, "WeakStrength") == 0){
        sscanf(entries[x].val, "%d", &s1);
        if (s1 > SMAXW) {
            printf("Maximum strength for weak items is %d!", SMAXW);
            s1 = SMAXW;
            }
           else if (s1 < 1) {
            printf("Minimum strength for weak items is 1!");
            s1 = 1;
            }
        }
    else if (strcmp(entries[x].name, "StrongStrength") == 0) {
        sscanf(entries[x].val, "%d", &s2);
        if (s2 > SMAXS) {
            printf("Maximum strength for strong items is %d!", SMAXS);
            s2 = SMAXS;
            }
        else if (s2 < 2) {
            printf("Minimum strength for strong items is 2!");
            s2 = 2;
            }
        }
    else if (strcmp(entries[x].name, "NumberOfFeatures") == 0) {
        sscanf(entries[x].val, "%d", &features);
        if (features > 32) {
            printf("Maximum number of features is 32!");
            features = 32;
            }
        else if (features < 4) {
            printf("Minimum number of features is 4!");
            features = 4;
            }
        }
    else if (strcmp(entries[x].name, "RecognitionCriterion") == 0)
        sscanf(entries[x].val, "%f", &crit);
    else if (strcmp(entries[x].name, "RecallCriterion") == 0)
        sscanf(entries[x].val, "%f", &recall_crit);
    else if (strcmp(entries[x].name, "ProbabilityOfStorage") == 0)
        sscanf(entries[x].val, "%f", &learn);
    else if (strcmp(entries[x].name, "NumberOfSimulations") == 0){
        sscanf(entries[x].val, "%d", &nsim);
        if (nsim > NSIM) {
            printf("Maximum number simulatoins is %d!", NSIM);
            nsim = NSIM;
            }
        }
    }

    printf("<H1>Query Results</H1>");
    if (s2 <= s1)
        printf("Strong items should be stronger than weak ones!<BR>\n");

    printf("<PRE>\n");
    ranset();                                 /* init the random no generator */
    simulate();
    screen();
    printf("</PRE>\n");
}



/*********************
 * FUNCTION SIMULATE *
 *********************
 *
 * main simulation routine
 *
 */
void simulate(void)
{
	int i;

				count_cw  = 0;                         /* initialize the counters */
	count_cs  = 0;
	count_hw  = 0;
	count_hs  = 0;
	count_f   = 0;
	count_fs  = 0;
	count     = 0;

	famtw     = 0.0;
	famtw2    = 0.0;
	famts     = 0.0;
	famts2    = 0.0;
	famd      = 0.0;
	famd2     = 0.0;
	famds     = 0.0;
	famds2    = 0.0;

	corw      = 0.0;
	cors      = 0.0;
  mycorw    = 0.0;
  mycors    = 0.0;
																					 /* initialize the storage arrays */
	for (i=0; i<NSIM; i++) {
		VARTW[i]  = 0.0;
		VARTS[i]  = 0.0;
		VARD[i]   = 0.0;
		VARDS[i]  = 0.0;
		DPRPW[i]  = 0.0;
		DPRPS[i]  = 0.0;
		DPRDW[i]  = 0.0;
		DPRDS[i]  = 0.0;
	}

	if (strength==1)
		length=s1*n;
	else if (strength==2)
		length=s2*n;
	else
		length=(s1*n/2)+(s2*n/2);
	for (i=1; i<=nsim; i++) {
		making_matrix();                  /* initialize the matrix */
		study();
		switch(cond) {
			case 1:
				recognition_single();
				break;
			case 2:
				recognition_ass_rearr();
				break;
			case 3:
				recognition_ass_mixed();
				break;
			case 4:
				recognition_ass_new();
				break;
			case 5:
				cued_recall();
				break;
		}
	}
}

/**************************
 * FUNCTION MAKING_MATRIX *
 **************************
 *
 * This function sets up a new lexicon
 *
 */
void making_matrix(void)
{
	int i,j,k;
	int check,temp,choose;



	/* creates the features of each item in the lexicon */

	for (i=0; i<(2*n); i++) {
L1: for (k=0; k<(2*features); k++)
			items[i][k]=gau((double) mu, (double)sigma/(double)features);


	/* makes sure that each pair of items item is unique */

		for (j=0; j<i; j++) {
			check=0;
			for (k=0; k<features; k++)
				if (items[i][k]!=items[j][k])
					check=1;
			if (check==0)
				goto L1;
			check=0;
			for (k=features; k<2*features; k++)
				if (items[i][k]!=items[j][k])
					check=1;
			if (check==0)
				goto L1;
		}
	}
}

/******************
 * FUNCTION STUDY *
 ******************
 *
 * creation of memory due to list presentation
 *
 */
void study(void)
{
int i, j, x,y;


	for (x=0; x<features; x++)
		for (y=0; y<features; y++)
			new_memory[x][y] = 0.0;

	if (strength==1) {   /* pure weak */
		for (i=0; i<n; i++) {  /* n is the length of the list (# pairs) */
			for (x=0;x<features;x++)
				for (y=0;y<features;y++)
					matrix1[x][y] = items[i][x]*items[i][y+features];

			for (j=0;j<s1;j++)
				for (x=0;x<features;x++)
					for (y=0;y<features;y++)
						if ((float) randw(1000) / (float) 1000.0 < learn)
							new_memory[x][y] += matrix1[x][y];
						else
							new_memory[x][y] += gau ((double)mu, (double)sigma/(double)features) * gau ((double)mu, (double)sigma/(double)features);

		}
	}
	else if (strength==2) {   /* pure strong */
		for (i=0; i<n; i++) {
			for (x=0;x<features;x++)
				for (y=0;y<features;y++)
					matrix1[x][y] = (float) items[i][x]*items[i][y+features];

			for (j=0;j<s2;j++)
				for (x=0;x<features;x++)
					for (y=0;y<features;y++)
						if ((float) randw(1000) / (float) 1000.0 < learn)
							new_memory[x][y] += matrix1[x][y];
						else
							new_memory[x][y] += gau ((double)mu, (double)sigma/(double)features) * gau ((double)mu, (double)sigma/(double)features);

		}
	}
	else {       /* mixed list */
		for (i=0; i<n/2; i++) {
			for (x=0;x<features;x++)
				for (y=0;y<features;y++)
					matrix1[x][y] = (float) items[i][x]*items[i][y+features];

			for (j=0;j<s1;j++)
				for (x=0;x<features;x++)
					for (y=0;y<features;y++)
						if ((float) randw(1000) / (float) 1000.0 < learn)
							new_memory[x][y] += matrix1[x][y];
						else
							new_memory[x][y] += gau ((double)mu, (double)sigma/(double)features) * gau ((double)mu, (double)sigma/(double)features);


		}
		for (i=n/2; i<n; i++) {
			for (x=0;x<features;x++)
				for (y=0;y<features;y++)
					matrix1[x][y] = (float) items[i][x]*items[i][y+features];

			for (j=0;j<s2;j++)
				for (x=0;x<features;x++)
					for (y=0;y<features;y++)
						if ((float) randw(1000) / (float) 1000.0 < learn)
							new_memory[x][y] += matrix1[x][y];
						else
							new_memory[x][y] += gau ((double)mu, (double)sigma/(double)features) * gau ((double)mu, (double)sigma/(double)features);

		}
	}

}

/*******************************
 * FUNCTION RECOGNITION_SINGLE *
 *******************************
 *
 * single item recognition simulation routine
 *
 */
void recognition_single(void)
{
	static real hrate, farate, dp;
	extern int dprime_();
	int i;
	float t;
	float fam;                        /* familiarity of test-item */
	int chw  = 0,
			chs  = 0,
			cf   = 0;
	float ftw   = 0.0,          /* variables to compute means and variances */
				ftw2  = 0.0,          /* for the target and distractor distibutions */
				fts   = 0.0,
				fts2  = 0.0,
				fd    = 0.0,
				fd2   = 0.0,
				vtw   = 0.0,
				vts   = 0.0,
				vdw   = 0.0;

	/* targets */

	for (i=0; i<n; i++) {
		/* checks for recognition of the item on the x axis */
		fam=single_fam1(i);
		/* checks for recognition of the item on the y axis */
		fam+=single_fam2(i);
		if ((strength == 1) || (strength == 3 && i < n/2)) {
			famtw += fam;
                        printf("%f\n", fam);
			famtw2 += fam * fam;
			ftw += fam;
			ftw2 += fam * fam;
		}
		if ((strength == 2) || (strength == 3 && i >= n/2)) {
			famts += fam;
			famts2 += fam * fam;
			fts += fam;
			fts2 += fam * fam;
		}
		if (fam > crit) {   /* familiarity above recognition criterion threshold*/
			if ((strength == 1) || (strength == 3 && i < n/2)) {
				count_hw++;
				chw++;
			}
			if ((strength == 2) || (strength == 3 && i >= n/2)) {
				count_hs++;
				chs++;
			}
		}
	}
	/* distractors */
	for (i=n; i<(2*n); i++) {
		/* checks for recognition of the item on the x axis */
		fam = single_fam1(i);
                printf("%f\n", fam);
		/* checks for recognition of the item on the y axis */
		fam += single_fam2(i);

		famd += fam;
		famd2 += fam * fam;
		fd += fam;
		fd2 += fam * fam;

		if (fam > crit) {
			count_f++;
			cf++;
		}
	}

	/* compute results and store in arrays */
	if (strength == 3)
		t = (float)n/2;
	else
		t = (float)n;


	vtw = ((t*ftw2)-(ftw*ftw))/(t*(t-1.0));
	if (vtw < 0.0001)
		vtw = 0.0001;
	VARTW[count]  = vtw;

	vts = ((t*fts2)-(fts*fts))/(t*(t-1.0));
	if (vts < 0.0001)
		vts = 0.0001;
	VARTS[count]  = vts;

	vdw = ((n*fd2)-(fd*fd))/(n*(n-1.0));
	if (vdw < 0.0001)
		vdw = 0.0001;
	VARD[count]   = vdw;

	hrate = (float)chw/t;
	farate = (float)cf/(float)n;
	dprime_(&hrate, &farate, &dp);
	DPRPW[count]  = dp;

	hrate = (float)chs/t;
	farate = (float)cf/(float)n;
	dprime_(&hrate, &farate, &dp);
	DPRPS[count]  = dp;

	DPRDW[count]  = ((ftw/t)-(fd/(float)n))/sqrt(VARD[count]);
	DPRDS[count]  = ((fts/t)-(fd/(float)n))/sqrt(VARD[count]);

	count++;
}

/**********************************
 * FUNCTION RECOGNITION_ASS_REARR *
 **********************************
 *
 * associative recognition simulation routine
 * targets    : A-B (intact)
 * distractors: A-D (rearranged)
 *
 */
void recognition_ass_rearr(void)
{
	static real hrate, farate, dp;
	extern int dprime_();
	int i, j;
	float t;
	float fam;                                    /* familiarity of test-item */
	int chw   = 0,
			chs   = 0,
			cf    = 0,
			cfs   = 0;
	float ftw   = 0.0,          /* variables to compute means and variances */
				ftw2  = 0.0,          /* for the target and distractor distibutions */
				fts   = 0.0,
				fts2  = 0.0,
				fd    = 0.0,
				fd2   = 0.0,
				fds   = 0.0,
				fds2  = 0.0,
				vtw   = 0.0,
				vts   = 0.0,
				vdw   = 0.0,
				vds   = 0.0;

	for (i=0; i<n/4; i++) {                                     /* targets 1 */
		fam = double_fam(i,i);
		if (strength == 2) {
			famts += fam;
			famts2 += fam * fam;
			fts += fam;
			fts2 += fam * fam;
		}
		else {
			famtw += fam;
			famtw2 += fam * fam;
			ftw += fam;
			ftw2 += fam * fam;
		}

		if (fam > crit) {
			if (strength == 2) {
				count_hs++;
				chs++;
			}
			else {
				count_hw++;
				chw++;
			}
		}
	}

	for (i=n/4; i<n/2; i++) {                   /* distractors 1 */
		j=i+1;
		if (i==((n/2)-1))
			j=n/4;
		fam = double_fam(i,j);
		famd += fam;
		famd2 += fam * fam;
		fd += fam;
		fd2 += fam * fam;

		if (fam > crit) {
			count_f++;
			cf++;
		}
	}

	for (i=n/2; i<3*n/4; i++) {                               /* targets 2 */
		fam = double_fam(i,i);
		if (strength == 1) {
			famtw += fam;
			famtw2 += fam * fam;
			ftw += fam;
			ftw2 += fam * fam;
		}
		else {
			famts += fam;
			famts2 += fam * fam;
			fts += fam;
			fts2 += fam * fam;
		}

		if (fam > crit) {
			if (strength == 1) {
				count_hw++;
				chw++;
			}
			else {
				count_hs++;
				chs++;
			}
		}
	}

	for (i=3*n/4; i<n; i++) {                 /* distractors 2 */
		j=i+1;
		if (i==n-1)
			j=3*n/4;
		fam = double_fam(i,j);
		if (strength == 3) {
			famds += fam;
			famds2 += fam * fam;
			fds += fam;
			fds2 += fam * fam;
		}
		else {
			famd += fam;
			famd2 += fam * fam;
			fd += fam;
			fd2 += fam * fam;
		}

		if (fam > crit) {
			if (strength == 3) {
				count_fs++;
				cfs++;
			}
			else {
				count_f++;
				cf++;
			}
		}
	}
																		 /* compute results and store in arrays */
	if (strength == 3)
		t = (float)n/4;
	else
		t = (float)n/2;

	vtw = ((t*ftw2)-(ftw*ftw))/(t*(t-1));
	if (vtw < 0.0001)
		vtw = 0.0001;
	VARTW[count]  = vtw;

	vts = ((t*fts2)-(fts*fts))/(t*(t-1));
	if (vts < 0.0001)
		vts = 0.0001;
	VARTS[count]  = vts;

	vdw = ((t*fd2)-(fd*fd))/(t*(t-1));
	if (vdw < 0.0001)
		vdw = 0.0001;
	VARD[count]   = vdw;

	if (strength == 3) {
		vds = ((t*fds2)-(fds*fds))/(t*(t-1));
		if (vds < 0.0001)
			vds = 0.0001;
		VARDS[count]  = vds;
	}

	hrate = (float)chw/t;
	farate = (float)cf/t;
	dprime_(&hrate, &farate, &dp);
	DPRPW[count]  = dp;

	hrate = (float)chs/t;
	if (strength == 3 )
		farate = (float)cfs/t;
	else
		farate = (float)cf/t;
	dprime_(&hrate, &farate, &dp);
	DPRPS[count]  = dp;

	DPRDW[count]  = ((ftw/t)-(fd/t))/sqrt(VARD[count]);
	if (strength == 3 )
		DPRDS[count]  = ((fts/t)-(fds/t))/sqrt(VARDS[count]);
	else
		DPRDS[count]  = ((fts/t)-(fd/t))/sqrt(VARD[count]);

	count++;
}

/**********************************
 * FUNCTION RECOGNITION_ASS_MIXED *
 **********************************
 *
 * associative recognition simulation routine
 * targets    : A-B (intact)
 * distractors: A-X (mixed)
 *
 */
void recognition_ass_mixed(void)
{
	static real hrate, farate, dp;
	extern int dprime_();
	int i, j;
	float t;
	float fam;                                    /* familiarity of test-item */
	int chw   = 0,
			chs   = 0,
			cf    = 0,
			cfs   = 0;
	float ftw   = 0.0,          /* variables to compute means and variances */
				ftw2  = 0.0,          /* for the target and distractor distibutions */
				fts   = 0.0,
				fts2  = 0.0,
				fd    = 0.0,
				fd2   = 0.0,
				fds   = 0.0,
				fds2  = 0.0,
				vtw   = 0.0,
				vts   = 0.0,
				vdw   = 0.0,
				vds   = 0.0;

	for (i=0; i<n/4; i++) {                                     /* targets 1 */
		fam = double_fam(i,i);
		if (strength == 2) {
			famts += fam;
			famts2 += fam * fam;
			fts += fam;
			fts2 += fam * fam;
		}
		else {
			famtw += fam;
			famtw2 += fam * fam;
			ftw += fam;
			ftw2 += fam * fam;
		}

		if (fam > crit) {
			if (strength == 2) {
				count_hs++;
				chs++;
			}
			else {
				count_hw++;
				chw++;
			}
		}
	}

	for (i=n/4; i<n/2; i++) {                   /* distractors 1 */
		j=n+i;
		fam = double_fam(i,j);
		famd += fam;
		famd2 += fam * fam;
		fd += fam;
		fd2 += fam * fam;

		if (fam > crit) {
			count_f++;
			cf++;
		}
	}

	for (i=n/2; i<3*n/4; i++) {                               /* targets 2 */
		fam = double_fam(i,i);
		if (strength == 1) {
			famtw += fam;
			famtw2 += fam * fam;
			ftw += fam;
			ftw2 += fam * fam;
		}
		else {
			famts += fam;
			famts2 += fam * fam;
			fts += fam;
			fts2 += fam * fam;
		}

		if (fam > crit) {
			if (strength == 1) {
				count_hw++;
				chw++;
			}
			else {
				count_hs++;
				chs++;
			}
		}
	}

	for (i=3*n/4; i<n; i++) {                 /* distractors 2 */
		j=(3*n/2)+i;
		fam = double_fam(i,j);
		if (strength == 3) {
			famds += fam;
			famds2 += fam * fam;
			fds += fam;
			fds2 += fam * fam;
		}
		else {
			famd += fam;
			famd2 += fam * fam;
			fd += fam;
			fd2 += fam * fam;
		}

		if (fam > crit) {
			if (strength == 3) {
				count_fs++;
				cfs++;
			}
			else {
				count_f++;
				cf++;
			}
		}
	}
																		 /* compute results and store in arrays */
	if (strength == 3)
		t = (float)n/4;
	else
		t = (float)n/2;

	vtw = ((t*ftw2)-(ftw*ftw))/(t*(t-1));
	if (vtw < 0.0001)
		vtw = 0.0001;
	VARTW[count]  = vtw;

	vts = ((t*fts2)-(fts*fts))/(t*(t-1));
	if (vts < 0.0001)
		vts = 0.0001;
	VARTS[count]  = vts;

	vdw = ((t*fd2)-(fd*fd))/(t*(t-1));
	if (vdw < 0.0001)
		vdw = 0.0001;
	VARD[count]   = vdw;

	if (strength == 3) {
		vds = ((t*fds2)-(fds*fds))/(t*(t-1));
		if (vds < 0.0001)
			vds = 0.0001;
		VARDS[count]  = vds;
	}

	hrate = (float)chw/t;
	farate = (float)cf/t;
	dprime_(&hrate, &farate, &dp);
	DPRPW[count]  = dp;

	hrate = (float)chs/t;
	if (strength == 3 )
		farate = (float)cfs/t;
	else
		farate = (float)cf/t;
	dprime_(&hrate, &farate, &dp);
	DPRPS[count]  = dp;

	DPRDW[count]  = ((ftw/t)-(fd/t))/sqrt(VARD[count]);
	if (strength == 3 )
		DPRDS[count]  = ((fts/t)-(fds/t))/sqrt(VARDS[count]);
	else
		DPRDS[count]  = ((fts/t)-(fd/t))/sqrt(VARD[count]);

	count++;
}

/********************************
 * FUNCTION RECOGNITION_ASS_NEW *
 ********************************
 *
 * associative recognition simulation routine
 * targets    : A-B
 * distractors: X-Y
 *
 */
void recognition_ass_new(void)
{
	static real hrate, farate, dp;
	extern int dprime_();
	int i, j;
	float t;
	float fam;                                    /* familiarity of test-item */
	int chw   = 0,
			chs   = 0,
			cf    = 0;
	float ftw   = 0.0,          /* variables to compute means and variances */
				ftw2  = 0.0,          /* for the target and distractor distibutions */
				fts   = 0.0,
				fts2  = 0.0,
				fd    = 0.0,
				fd2   = 0.0,
				vtw   = 0.0,
				vts   = 0.0,
				vdw   = 0.0;

	for (i=0; i<n/2; i++) {                                     /* targets 1 */
		fam = double_fam(i,i);
		if (strength == 2) {
			famts += fam;
			famts2 += fam * fam;
			fts += fam;
			fts2 += fam * fam;
		}
		else {
			famtw += fam;
			famtw2 += fam * fam;
			ftw += fam;
			ftw2 += fam * fam;
		}

		if (fam > crit) {
			if (strength == 2) {
				count_hs++;
				chs++;
			}
			else {
				count_hw++;
				chw++;
			}
		}
	}

	for (i=n; i<2*n; i++) {                   /* distractors */
		fam = double_fam(i,i);
		famd += fam;
		famd2 += fam * fam;
		fd += fam;
		fd2 += fam * fam;

		if (fam > crit) {
			count_f++;
			cf++;
		}
	}

	for (i=n/2; i<n; i++) {                               /* targets 2 */
		fam = double_fam(i,i);
		if (strength == 1) {
			famtw += fam;
			famtw2 += fam * fam;
			ftw += fam;
			ftw2 += fam * fam;
		}
		else {
			famts += fam;
			famts2 += fam * fam;
			fts += fam;
			fts2 += fam * fam;
		}

		if (fam > crit) {
			if (strength == 1) {
				count_hw++;
				chw++;
			}
			else {
				count_hs++;
				chs++;
			}
		}
	}
	if (strength == 3)
		t = (float)n/2;
	else
		t = (float)n;

	vtw = ((t*ftw2)-(ftw*ftw))/(t*(t-1));
	if (vtw < 0.0001)
		vtw = 0.0001;
	VARTW[count]  = vtw;

	vts = ((t*fts2)-(fts*fts))/(t*(t-1));
	if (vts < 0.0001)
		vts = 0.0001;
	VARTS[count]  = vts;
	vdw = (((float)n * fd2)-(fd*fd))/((float)n*(float)n-1);
	if (vdw < 0.0001)
		vdw = 0.0001;
	VARD[count]   = vdw;

	hrate = (float)chw/t;
	farate = (float)cf/(float)n;
	dprime_(&hrate, &farate, &dp);
	DPRPW[count]  = dp;

	hrate = (float)chs/t;
	farate = (float)cf/(float)n;
	dprime_(&hrate, &farate, &dp);
	DPRPS[count]  = dp;

	DPRDW[count]  = ((ftw/t)-(fd/(float)n))/sqrt(VARD[count]);
	DPRDS[count]  = ((fts/t)-(fd/(float)n))/sqrt(VARD[count]);

	count++;
}



/************************
 * FUNCTION CUED_RECALL *
 ************************
 *
 * cued recall simulation routine
 * A-B1..n & A1..n-B
 *
 * The target A is paired with every second item in the target pairs and
 * the dot product is taken of each combination.  The combination with the
 * highest dot_product is retained, and if it is the correct target pair,
 * correct recall has occured.
 *
 */
void cued_recall (void)
{
int   j,k,x,y, recall_item;
float	fam, max_fam,
			total;



	for (j=0;j<n;j++) {
		max_fam=0.0;
    total=0.0;
		for (k=0;k<n;k++) {
			/* make an association matrix of j * k */

			fam = double_fam (j, k);

			total += fam;

			if (fam > max_fam) {
				max_fam = fam;
				recall_item = k;
			}
		}

		if (strength==1) {
			corw+=total/((float)n);
		}
		else
			if (strength==2) {
				cors+=total/((float)n);
			}
			else {
				if (j<n/2) {
					corw+=total/((float)n);
				}
				else {
					cors+=total/((float)n);
				}
			}

		if (recall_item == j && max_fam > recall_crit) {
			if (strength == 1)
				count_cw++;
			if (strength == 2)
				count_cs++;
			if (strength == 3 && j<n/2)

				count_cw++;
			if (strength == 3 && j>=n/2)
				count_cs++;
		}
	}
}



/*******************
 * FUNCTION RANSET *
 *******************
 *
 * set up the random seed
 */
void ranset(void)
{
	struct timeb timebuffer;

	ftime(&timebuffer);
	srand(timebuffer.millitm);

/*	randomize ();  *//* this can be used in turbo c instead */
}


/*****************
 * FUNCTION RAND *
 *****************
 *
 * random number generator between 0 and 1
 *
 */
float rand1(void)
{
	return((float)rand()/RAND_MAX);
}

/****************
 * FUNCTION GAU *
 ****************
 *
 * gaussian distribution, truncated at 0.0
 *
 */
float gau(double mu, double sigma)
{
	float random_array[12];
	int i;
	float total, result;


	total =0.0;
	for (i=0; i<=11; i++) {
		random_array[i] = rand1();
		total += random_array[i];
	}
	result = (sqrt(sigma) * ( total - 6.0)) + mu;


	return(result);

}

/*****************
 * function PERM *
 *****************
 *
 * Takes and array (list) of length (len) and permutes the order of its
 * elements.
 *
 * The length (len) from PERM() becomes the range delimiter (n) in RANDW()
 */
void perm(int *list, int len)
{
	int i, j, temp;

	for (i = 0; i < len; i++) {
		j = randw(len);
		temp = list[j];
		list[j] = list[i];
		list[i] = temp;
	}

} /* end function PERM */

/******************
 * function RANDW *
 ******************
 *
 * Calls the random number generator RAND() and uses the modulus (%)
 * operator to limit the random numbers to a specific range
 */
int randw(int n)
{
	return(rand()%n);

} /* end function RANDW */

/*******************
 * FUNCTION DPRIME *
 *******************
 *
 * computes d-prime based on hit and false alarm rate
 *
 */
int dprime_(hr, far1, d)

real *hr, *far1, *d;

{
		static real y, z1, z2, y2;
		extern /* Subroutine */ int znorma_();

	if (*hr == (float)1) {
		*hr = (float).99;
	}
	if (*hr == (float)0) {
		*hr = (float).01;
	}
	if (*far1 == (float)1) {
		*far1 = (float).99;
	}
	if (*far1 == (float)0) {
		*far1 = (float).01;
	}
	znorma_(hr, &z1, &y);
	znorma_(far1, &z2, &y2);
	*d = z1 - z2;

	return 0;

} /* dprime_ */

/* Subroutine */
int znorma_(p, zz, y)
real *p, *zz, *y;
{
		/* Initialized data */

		static real c1 = (float)2.515517;
		static real c2 = (float).802853;
		static real c3 = (float).010328;
		static real c4 = (float)1.432788;
		static real c5 = (float).189269;
		static real c6 = (float).001308;

		/* System generated locals */
		real r_1;

		/* Builtin functions */
/*    double log(), sqrt(), r_sign(), exp();*/

		/* Local variables */
		static real t, ab, cc, dd, pq, sz;

		if ((r_1 = *p - (float).5) < (float)0.) {
	goto L11;
		} else if (r_1 == 0.0) {
	goto L10;
		} else {
	goto L9;
		}
L9:
		pq = (float)1. - *p;
		sz = (float)1.;
		goto L12;
L10:
		*zz = (float)0.;
		goto L99;
L11:
		pq = *p;
		sz = (float)-1.;
L12:
		t = sqrt(log(pq) * (-2.0));
		ab = ((c6 * t + c5) * t + c4) * t + (float)1.;
		cc = (c3 * t + c2) * t + c1;
		dd = cc / ab;
		*zz = t - dd;
		*zz = r_sign(zz, &sz);
L99:
		*y = exp(-(double)(*zz) * (*zz) / (float)2.) * (float).39894228;
		return 0;
} /* znorma_ */



double r_sign(a, b)
real *a, *b;
{
	if (*b<0.0) return(-1.0*fabs(*a));
	else return(fabs(*a));
}

/********************
 * FUNCTION ACTION2 *
 ********************
 *
 * performs action for MENU2 based on cursor position
 * Type of test
 */
void action2(int pos)
{
	switch(pos) {
		case 0:
			cond = 1;
			calculate_criterion ();
			break;
		case 1:
			cond = 2;
			calculate_criterion ();
			break;
		case 2:
			cond = 3;
			calculate_criterion ();
			break;
		case 3:
			cond = 4;
			calculate_criterion ();
			break;
		case 4:
			cond = 5;
			calculate_criterion ();
			break;
		case 5:
			break;
	}
}

/*********************
 * FUNCTION SETTINGS *
 *********************
 *
 * displays values for variables under the menus
 *
 */
void settings(void)
{
	printf("Current settings:\n");
	printf("\tTest:\t\t\t\t%s (%d sim.)\n", menu2[cond-1], nsim);
	printf("\tLength:\t\t\t\t%d\n", n);
	if (strength == 1)
		printf("\tStrength:\t\t\t%s; value = %d\n", menu4[strength-1], s1);
	if (strength == 2)
		printf("\tStrength:\t\t\t%s; value = %d\n", menu4[strength-1], s2);
	if (strength == 3)
		printf("\tStrength:\t\t\t%s; weak = %d, strong = %d\n",
							menu4[strength-1], s1, s2);

	printf("\t%s:\t\tcrit = %3.3f\n", menu5[0], crit);
	printf("\t%s:\t\tcrit = %3.3f \n", menu5[1], recall_crit);
	printf("\t%s:\t\tfeatures = %d\n", menu5[2], features);
	printf("\t%s:\t\tlearn (L) = %.2f\n", menu5[3], learn);
	printf("\t%s:\t\t\tmean = %3.2f\n", menu5[4], mu);
	printf("\t%s:\t\t\tpower = %3.2f\n", menu5[5], sigma);

}

/*******************
 * FUNCTION SCREEN *
 *******************
 *
 * put settings() and results of simulation on screen
 *
 */
void screen(void)
{
	char junk='\0';
	int i;
	float t, N, S;
	static real hrate, farate, dp;
	extern int dprime_();

																								 /* intialize the variables */
	pcw = pcs = 0.0;
	hw = mtw = vartw = tvartw = 0.0;
	hs = mts = varts = tvarts = 0.0;
	f  = md  = vard  = tvard  = 0.0;
	fs = mds = vards = tvards = 0.0;
	dprpw = tdprpw = dprps = tdprps = 0.0;
	dprdw = tdprdw = dprds = tdprds = 0.0;


	settings();
	N = (float)nsim;


	if (cond == 2 || cond == 3) {
		if (strength == 3)
			t = (float)n/4.0;
		else
			t = (float)n/2.0;
	}
	else {
		if (strength == 3)
			t = (float)n/2.0;
		else
			t = (float)n;
	}

	S = t*N;


	if (cond == 5) {                                   /* results cued recall */
		pcw=(float)count_cw/(N*t);
		pcs=(float)count_cs/(N*t);
		if (strength == 3) {
			printf("\n\nWeak items:\tP(C) = %.3f\n", pcw);
			printf ("Average correlation = %1.3f\n", mycorw/(float)nsim);
			printf("Mean familiarity = %5.2f\n\n",corw/((float)nsim*(float)n));
      printf("Strong items:\tP(C) = %.3f\n", pcs);
			printf ("Average correlation = %1.3f\n", cors/S);
		}
		else if (strength == 2) {
			printf("\n\nStrong items:\tP(C) = %.3f\n", pcs);
			printf ("Average correlation = %1.3f\n", cors/S);
		}
		else {
			printf("\n\nWeak items:\tP(C) = %.3f\n", pcw);
			printf ("Average correlation = %1.3f\n", corw/S);
		}
	}

	else {                                             /* results recognition */
		for (i=0; i<N; i++) {
			vartw += VARTW[i];
			varts += VARTS[i];
			vard  += VARD[i];
			vards += VARDS[i];
			dprpw += DPRPW[i];
			dprps += DPRPS[i];
			dprdw += DPRDW[i];
			dprds += DPRDS[i];
		}
		vartw   /= N;
		varts   /= N;
		vard    /= N;
		vards   /= N;
		dprpw   /= N;
		dprps   /= N;
		dprdw   /= N;
		dprds   /= N;

		hw      = (float)count_hw/S;
		mtw     = famtw/S;
		tvartw  = ((S*famtw2)-(famtw*famtw))/(S*(S-1));

		hs      = (float)count_hs/S;
		mts     = famts/S;
		tvarts  = ((S*famts2)-(famts*famts))/(S*(S-1));

		if (cond == 1 || cond == 4) {
			t     = (float)n;
			S     = t*N;
		}

		f       = (float)count_f/S;
		md      = famd/S;

		tvard   = ((S*famd2)-(famd*famd))/(S*(S-1.0));

		if ((cond == 2 || cond == 3) && strength == 3) {
			fs     = (float)count_fs/S;
			mds    = famds/S;

			tvards = ((S*famds2)-(famds*famds))/(S*(S-1));
		}

		hrate   = hw;
		farate  = f;

		dprime_(&hrate, &farate, &dp);
		tdprpw  = dp;

		hrate   = hs;
		if ((cond == 2 || cond == 3) && strength == 3)
			farate  = fs;
		else
			farate  = f;
		dprime_(&hrate, &farate, &dp);
		tdprps  = dp;
		if (tvard < 0.0001)
			tvard = 0.0001;
		if (tvards < 0.001)
			tvards = 0.0001;
		tdprdw  = (mtw-md)/sqrt(tvard);
		if (tdprdw < 0.0001)
			tdprdw = 0.0001;
		if ((cond == 2 || cond == 3) && strength == 3)
			tdprds  = (mts-mds)/sqrt(tvards);
		else
			tdprds  = (mts-md)/sqrt(tvard);

		printf("\n\nResults averaged over simulations, and for total simulation"
						 " (in brackets)\n");
		if (strength == 3)
			printf ("\t\tTargets\t\t\t\t\t\tDistractors\n");
		else
			printf("\t\tTargets\t\t\t\tDistractors\n");
		if (strength == 1) {
			printf("\nproportion\t%10.3f\t\t\t%10.3f\n", hw,f);
			printf("mean fam.\t%10.3f\t\t\t%10.3f\n", mtw, md);
			printf("variance \t%10.3f (%10.3f)\t\t%10.3f (%10.3f)\n\n",
						vartw, tvartw, vard, tvard);
			printf("d' (prop.)\t%.3f (%.3f)\n\n", dprpw, tdprpw);
			printf("d' (distr.)\t%.3f (%.3f)\n", dprdw, tdprdw);
		}
		if (strength == 2) {
			printf("\nproportion\t%10.3f\t\t\t%10.3f\n", hs,f);
			printf("mean fam.\t%10.3f\t\t\t%10.3f\n", mts, md);
			printf("variance \t%10.3f (%10.3f)\t\t%10.3f (%10.3f)\n\n",
						varts, tvarts, vard, tvard);
			printf("d' (prop.)\t%.3f (%.3f)\n\n", dprps, tdprps);
			printf("d' (distr.)\t%.3f (%.3f)\n", dprds, tdprds);
		}
		if (strength == 3) {
			if (cond == 2 || cond == 3)
				printf("\t\tweak\t\tstrong\t\tweak\t\tstrong\n");
			else
				printf("\t\tweak\t\t\t  strong\n");

			if (cond == 2 || cond == 3) {
				printf("proportion%7.2f\t\t%7.2f\t\t%7.2f\t\t%7.2f\n", hw, hs, f, fs);
				printf("mean fam  %7.2f\t\t%7.2f\t\t%7.2f\t\t%7.2f\n", mtw, mts, md, mds);
				printf("variance  %7.2f(%7.2f)\t%7.2f(%7.2f)%7.2f(%7.2f)%7.2f(%7.2f)"
							"\n\n",vartw, tvartw,varts,tvarts,vard,tvard,vards,tvards);
			}
			else {
				printf("proportion\t%9.3f\t\t%9.3f\t\t%9.3f\n", hw, hs, f);
				printf("mean fam.\t%9.3f\t\t%9.3f\t\t%9.3f\n", mtw, mts, md);
				printf("variance \t%9.3f (%9.3f)\t%9.3f (%9.3f)\t%9.3f (%9.3f)\n\n",
							vartw, tvartw, varts, tvarts, vard, tvard);
			}

			printf("d' (prop.)\t%.3f (%.3f)\t%.3f (%.3f)\n\n",
						dprpw, tdprpw, dprps, tdprps);
			printf("d' (distr.)\t%.3f (%.3f)\t%.3f (%.3f)\n",
						dprdw, tdprdw, dprds, tdprds);
		}
	}

}


/***********************************************************************
*  Name:    single_fam1
*  params:  int i - an item in the first position of a pair (eg A- )
*  returns: float fam - the result of the dot product
*  desc:    This function creates an assocation matrix with an
*  item from the first position of a pair and another vector of all 1s.
*  eg  i vector is 1 0 1 1 so
*            _
*            1             1111
*            0             0000
*            1 * |1111|  = 1111
*            1             1111
*            -
* then the dot product of the resultant matrix with the memory is taken.
*
**/
float single_fam1 (int i)
{

int   x,y;
float fam=0.0;


	/* make a matrix of item*1 */
	for (x=0;x<features;x++)
		for (y=0;y<features;y++)
			matrix1[x][y] = items[i][x];


	fam = dot_product (matrix1, new_memory, features);


	return fam;

}


/***********************************************************************
*  Name:    single_fam2
*  params:  int i - an item in the second position of a pair (eg  -B)
*  returns: float fam - the result of the dot product
*  desc:    This function creates an assocation matrix with an
*  item from the second position of a pair and another vector of all 1s.
*  eg  i vector is 1 0 1 1 so
*            _
*            1             1011
*            1             1011
*            1 * |1011|  = 1011
*            1             1011
*            -
* then the dot product of the resultant matrix with the memory is taken.
*
**/
float single_fam2 (int i)
{

int   x,y;
float fam=0.0;

	/* make a matrix of 1*item */
	for (x=0;x<features;x++)
		for (y=0;y<features;y++)
			matrix1[x][y] = items[i][y+features];


	fam = dot_product (matrix1, new_memory, features);


	return fam;

}




/***********************
 * FUNCTION DOUBLE_FAM *
 ***********************
 *
 * for a double item probe, a familiarity is computed.
 *
 * A 2D matrix is created using an item in the first pair
 * position, combined with an item in the second pair
 * position.
 *
 */
float double_fam(int first, int second)
{
int   x,y;
float fam=0.0;

	/* make an association matrix of first * second */
	for (x=0;x<features;x++)
		for (y=0;y<features;y++)
			matrix1[x][y] = items[first][x]*items[second][y+features];

	fam = dot_product (matrix1, new_memory, features);

	return fam;

}



/***********************
 * FUNCTION DOT_PRODUCT*
 ***********************
 *
 * Computes the dot product of two 2D arrays of size columns and rows.
 *
 */
float dot_product (float matrix1[32][32],
									 float matrix2[32][32], int size)
{
int x, y;
double sum=0.0;


	for (x=0;x<size;x++)
		for (y=0;y<size;y++)
			sum += (double)matrix1[x][y]*(double)matrix2[x][y];

	return ((float)sum);

}



void calculate_criterion ()
{
float target, distr;
float m2, var, f2,u4,u3;


	f2  = features*features;
	var = sigma/features;
	m2  = var + mu*mu;
	u4  = (float) pow ((double) mu, 4.0);
	u3  = (float) pow ((double) mu, 3.0);

	switch (cond) {
		case 1:
			target = f2 * (learn*m2*mu + (1.0-learn)*u4) +
							 f2*(learn*u3 + (1.0-learn)*u4) +
							 (n-1)*2*(f2*(learn*u3 + (1.0-learn)*u4));
			distr  = n*(f2 * (learn*u3 + (1.0-learn)*u4) +
							 f2*(learn*u3 + (1.0-learn)*u4));
/*			distr  = n*2.0*(f2*(learn*u3+(1.0-learn)*u4));*/
			crit   = (target+distr)/2.0;
			break;
		case 2:
			target = f2*(learn*(float)pow((double)m2, 2.0) +
							 (1.0-learn)*u4) + (n-1)*f2*u4;

			distr  = f2*2.0*(learn*m2*mu*mu + (1.0-learn)*u4) + u4*(n-2)*f2;
			crit   = (target+distr)/2.0;
			break;
		case 3:
			target = f2*(learn*(float)pow((double)m2, 2.0) +
							 (1.0-learn)*u4) + (n-1)*f2*u4;

			distr  = f2*(learn*m2*mu*mu + (1.0-learn)*u4) +u4*(n-1)*f2;

			crit   = (target+distr)/2.0;
			break;
		case 4:
			target = f2*(learn*(float)pow((double)m2, 2.0) +
							 (1.0-learn)*u4) + (n-1)*f2*u4;

			distr  = f2*n*u4;
			crit   = (target+distr)/2.0;
			break;
		case 5:
			recall_crit = (f2*(learn*(float)pow((double)m2,2.0) +
										(1.0-learn)*u4) +
										(n-1)*f2*u4)/2;
			break;

	}
}