% RunExperimentCT
%
% Measure detection threshold for grating in noise
% Each call to RunExperimentCT will measure one threshold.
% This is meant to be called as a subprogram after parameters
% have been set, not to be executed directly.
%
% BUG: Student Matlab sometimes hangs up in the call to CmdWinToUpperLeft(),
% and continues if you hit Command-Period.

% 6/26/96  dgp  Wrote it.
% 3/8/97   dhb  Massive simplifications.
%               Much more documentation.
% 1/12/98  dgp	Noted problem in Student Matlab.
% 7/9/98   dgp  Cope when GetClut is not available.
% 4/08/02  awi  -Added platform conditionals to avoid SCREENWinToFront calls in Windows
%               -Added platform condtionals to use Matlab's Sound on Windows.
%               -Fixed GetClicks.dll so that we can count more than double clicks. 
% 4/13/02 dgp  Replaced platform-dependent call to obsolete SndPlay by unconditional call to Snd.
% 9/25/03 dhb  Force screen open in 32 bit mode.

% Open up the experimental screen
backgroundEntry = 128;
[window,screenRect]=Screen(screenNumber,'OpenWindow',backgroundEntry,[],32);
matlabClut = Screen(window,'GetClut'); % a misnomer. Screen 'OpenWindow' just set up the CLUT.
if isempty(matlabClut)
	matlabClut=ClutDefault(window);
end

if strcmp('PCWIN',computer)
    Screen(window,'Close');  
else
	% Bring command window to front and move it to upper left.
	% That way we can print lots of stuff in the command window
	% that tells the user what's going on.
	[cmdDx,cmdDy,bufferWindow] = CmdWinToUpperLeft(window,matlabClut);
end

% Provide introductory remarks.
fprintf(1,'\n**************************************\n\n');
fprintf('This program demonstrates how to measure a threshold.\n');
fprintf('In this case we are measuring the contrast required to detect a grating in noise.\n');
fprintf('After some initialization, we conduct a series of\n');
fprintf('two interval forced choice trials.  On each trial, your task is\n');
fprintf('to detect the grating.  The task is made difficult by the presence\n');
fprintf('of visual noise (%g RMS contrast).  If the grating was in the first\n',noiseRmsContrast);
fprintf('interval, click the mouse once.  If it was in the second interval,\n');
fprintf('click the mouse twice. (Click three times to quit.) The Quest \n');
fprintf('threshold seeking algorithm controls the contrast displayed on each \n');
fprintf('trial, to home in on your threshold.\n\n');
fprintf('We''ll now run %g trials.\n\n',nTrials);
fprintf('Click the mouse now to proceed\n');
GetClicks;
fprintf('\n');

% Seed random number generator
ClockRandSeed;

% Read screen calibration and compute cluts.
InitializeClutsCT;

% Auditory waveforms
beep=sin(2*pi*0.012*[0:900]);

% Make the signal waveforms
disp('	Making signals.');
MakeSignalsCT;

% Initialize Quest.  The threshold guess is set
% relative to the signal to noise ratio.
% We aim a bit high because observers do better if they
% clearly see the signal on the first few trials.
disp('	Initialize Quest parameters.');
nIntervals=2; nResponses = 2;
beta=3.5; delta=0.01; gamma=1/nIntervals;
pCorrect=0.82;
if (Enoise > 0)
	thresholdGuess = 0.5*log10(80/(Esignal/Enoise))-0.5;
else
	thresholdGuess = -1.5;
end
thresholdGuess = max(thresholdGuess,-1.5);
priorSd = 4;
q=QuestCreate(thresholdGuess,priorSd,pCorrect,beta,delta,gamma);

if strcmp('PCWIN',computer)
    [window,screenRect]=Screen(screenNumber,'OpenWindow',backgroundEntry);
else
	% Bring up screen window to show images.
	SCREENWinToFront(window,bufferWindow,cmdDx,cmdDy,offClut);
end

% Set the linear clut for the experiment
Screen(window,'SetClut',linearClut);
putRect = CenterRect(signalRect,screenRect);

% Set up text style for user interface.
Screen(window,'TextFont','Arial');
Screen(window,'TextSize',20);
Screen(window,'TextStyle',0);
Screen(window,'TextMode','srcOr');

% We should modify the CLUT to reserve entries 0 and 255 to show white and black, respectively,
% as required by the Mac OS 9.1+ Palette Manager. However, this hasn't been done yet.
black=255; % actually this produces white, but 0 produces gray.
white=255;

% Run the trials
theData = zeros(nTrials,4);
wrongRight = ['Wrong' ; 'Right'];

HideCursor;

for trial=1:nTrials
	
	% Make noise images
	MakeNoisesCT;
	
	% Get level to test
	trialCon=10^QuestQuantile(q); 
	trialCon=min(maxQuestCon,trialCon);

	% On first trial, wait for user.
	if trial==1
		fixationRect=CenterRect(SetRect(0,0,4,4),screenRect);
		Screen(window,'FillRect',white,fixationRect);
		instString = 'Click mouse when ready';
		Screen(window,'DrawText',instString,50,100,black);
		GetClicks;
		Screen(window,'DrawText',instString,50,100,backgroundEntry);
		instString = 'Click once for interval 1, twice for interval 2, three times to quit.';
		Screen(window,'DrawText',instString,50,100,black);
	end
	
	whichInterval=randi(nIntervals);
	signal0=zeros(size(signal1));
	for interval=1:nIntervals
		% Create stimulus by adding signal and noise
		noiseName=sprintf('noise%d',interval);
		if interval == whichInterval
			eval(['stimulus = trialCon*signal1+' noiseName ';']);
		else
			eval(['stimulus = trialCon*signal0+' noiseName ';']);
		end
		stimulus=backgroundEntry+entriesPerContrast*stimulus;

		% Display stimulus
		WaitSecs(0.3);
		Snd('Play',beep,22000);
		Screen(window,'PutImage',stimulus,putRect)
		Screen(window,'SetClut',linearClut);
		WaitSecs(durationMs/1000);
		Screen(window,'FillRect',backgroundEntry,putRect);
	end

	% Process response
	choice=GetClicks;
	if choice > 2
		break;
	end
	response=(choice==whichInterval);
	GiveFeedback(response);
	respString = sprintf('Trial %g: %s',trial,wrongRight(response+1,:));
	xLoc = screenRect(rectRight)/2-Screen(window,'TextWidth',respString)/2;
	yLoc = screenRect(rectBottom)/2;
	Screen(window,'DrawText',respString,xLoc,yLoc,black);
	WaitSecs(0.5);
	Screen(window,'DrawText',respString,xLoc,yLoc,backgroundEntry);
	WaitSecs(0.2);
	Screen(window,'FillRect',white,fixationRect);
	
	% Update Quest and store data
	q = QuestUpdate(q,log10(trialCon),response);
	theData(trial,1) = log10(trialCon);
	theData(trial,2) = response;
	theData(trial,3) = QuestMean(q);
	theData(trial,4) = QuestSd(q);
end

% Close up
ShowCursor;
Screen(window,'Close');

% Tell user about threshold
fprintf('\n**************************************\n\n');
fprintf('Threshold estimate is %g ± %g (mean ± sd, in log units).\n',QuestMean(q),QuestSd(q));
if 0 % boring
	fprintf('Array theData contains the data.  It has four\n');
	fprintf('columns.  The first column provides for each trial\n');
	fprintf('the log contrast value tested.  The second column tells\n');
	fprintf('whether the trial was right (1) or wrong (0)\n');
	fprintf('The third column gives Quest''s running threshold estimate.\n');
	fprintf('The fourth column gives Quest''s running estimate standard\n');
	fprintf('deviation.\n\n');
end
fprintf('The plot shows how the threshold estimate converged\n');
fprintf('over %g trials.\n\n',trial);

% Plot the threshold estimates.
plot(theData(1:trial,3),'g');
xlabel('Trial number');
ylabel('Log contrast threshold');
title('Threshold estimate vs. Trial Number');
drawnow;


