' Version 0.4   Date 4 Jan 2006
DECLARE SUB calcdr (outcome%(), order%, disc!, corr!)
DECLARE SUB setdwts ()
DECLARE FUNCTION randnorm! ()
DECLARE SUB writesumdatafile ()
DECLARE SUB writeeventdatafile (raceno%, outcome%(), eventno%)
DECLARE SUB calcdiscards (raceno%)
DECLARE SUB calcdiscardn (nr%)
DECLARE SUB showoverallstats ()
DECLARE SUB setability ()
DECLARE SUB eventstats (eventno%, raceno%)
DECLARE SUB raceres3 (raceno%)
DECLARE SUB raceres2 ()
DECLARE SUB raceres1 ()
DECLARE SUB standings (raceno%)
DECLARE SUB sort (n%, jibno%(), outcome%(), place!())
DECLARE SUB heatres (heatno%, n%, jibno%(), outcome%())
DECLARE SUB allocseed ()
DECLARE SUB makeheatsizes ()
DECLARE SUB initialparams ()

' The program simulates an HMS 2002 event.
' Boats are given an 'ability' between 0 and 1.
' There is a possibility of incidents in each heat, costing places.
' The number of boats to promote is a variable.
' Note that the jib number given to a boat is a proxy for its
'  expected final place in the standings.  For example, jib no "1"
'  is always the most able boat and should come first in any event.
' For speed and convenience, the measures of discrepancy between a
'  boat's actual place and her expected place are calculated by
'  comparing the boat's actual place against her jib number.
' Terminology
'  Event = a number of races
'  Race = a number of heats
'  Heat = a group of boats
'  In this simulation, heats are numbered 1, 2, ...
'  The lowest heat always runs first, so the winning boats can be
'   promoted to compete in the next heat

COMMON SHARED maxfleetsize%, minfleetsize%, fleetsize%
COMMON SHARED maxnumpromo%, numpromo%, maxnumraces%, numraces%
COMMON SHARED maxnumevents%, numevents%
COMMON SHARED debug%, progress%, showlist%
COMMON SHARED maxnumheats%, numheats%
COMMON SHARED maxheatsize%, minheatsize%, largestheatsize%
COMMON SHARED heatsize1%(), heatsize2%(), heatsize3%()
COMMON SHARED fleetjib%(), newfleet%()
COMMON SHARED abmethod%, abmethod$, ablist(), posnsd
COMMON SHARED lowerz, upperz, lowerzwt, middlezwt, upperzwt
COMMON SHARED incidplaceslost, possincid
COMMON SHARED pointslist&(), posnsum#(), posnsum2#()
COMMON SHARED correlx, correlx2, correlxy, wtpower, dpower, discrep, discreplimit
COMMON SHARED sumr, sumr2, sumd, sumd2
COMMON SHARED resarray%(), rarrayc(), darrayc(), rarrayi(), darrayi()
COMMON SHARED evarray()
COMMON SHARED discrate%, numdisc%
COMMON SHARED makefile%, makefilesum%, file$
COMMON SHARED dwts(), dwtmethod%, dwts$

ON ERROR GOTO handler

' Change these maxima to suit
maxfleetsize% = 80   ' Max number of boats
minfleetsize% = 12   ' Minimum number of boats
maxnumheats% = 5     ' Max number of heats
maxheatsize% = 20    ' Max number of boats per heat
minheatsize% = 8     ' Minimum number of boats per heat
maxnumraces% = 200   ' Max number of races to simulate
maxnumevents% = 200  ' Max number of events summary data to collect

DIM fleetjib%(maxfleetsize%)    ' Like a fleet board, shows the jib nos
				'  of the boats in the current race
DIM newfleet%(maxfleetsize%)    ' Used in allocating seedings for race 1
DIM ablist(maxfleetsize%)       ' Boat ability score in jib no order
DIM heatsize1%(maxnumheats%)    ' Number of boats in each race 1 heat
DIM heatsize2%(maxnumheats%)    ' Number of boats in each race 2 heat
DIM heatsize3%(maxnumheats%)    ' Number of boats in race 3+ heats
DIM pointslist&(maxfleetsize%)  ' Accumulative points tally of boats
DIM posnsum#(maxfleetsize%)     ' Becomes the average place of a boat
DIM posnsum2#(maxfleetsize%)    ' Becomes the standard deviation of boat places
DIM resarray%(maxfleetsize%, maxnumraces%)  ' The array of race-by-race results
DIM rarrayc(maxnumraces%)       ' Cumulative correlation r and discrep d
DIM darrayc(maxnumraces%)       '  values as event progresses -- these
				'  scores improve as the event goes on
DIM rarrayi(maxnumraces%)       ' Individual correlation r and discrep d
DIM darrayi(maxnumraces%)       '  values for a given race
DIM evarray(maxnumevents%, 4)   ' Array of r and d scores for the events
				'  in a simulation run
				'  1 Corr r, no discards
				'  2 Corr r after discards
				'  3 Discrep d, no discards
				'  4 Discrep d after discards
DIM dwts(maxfleetsize%)         ' Discrepancy weights

CLS

PRINT "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"
PRINT "                         REGATTA  SIMULATOR                      v0.4"
PRINT "                       (c)2005 Lester Gilbert"
PRINT "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"
PRINT
PRINT "The program simulates a number of races so the effects of changing"
PRINT "a variety of parameters can be studied, such as the number of boats"
PRINT "promoted and relegated in each heat, the number of discards allowed,"
PRINT "the maximum size of any heat, and so on.  Simulator data can be saved"
PRINT "to CSV files which import into Excel easily."
PRINT
PRINT "Current setup ";
PRINT "(sail numbers from 1 to"; maxfleetsize%; ")"
PRINT "is for up to"; maxnumraces%; "races per event with race by race results or"
PRINT "any number of races with recording of the first"; maxnumraces%; "results."
PRINT "Each race has no more than"; maxnumheats%; "heats and"; maxheatsize%; "yachts per heat."
PRINT "Up to"; maxnumevents%; "events can be repeated and summary data recorded,"
PRINT "or more than this with summary data for the first"; maxnumevents%; "recorded."
PRINT
PRINT "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++"

' Ask user for the simulator parameters
CALL initialparams
' CALL outputfileparams

' Construct the heat size arrays
CALL makeheatsizes

' Calculate some preparatory statistics
correlx = 0: correlx2 = 0
sumr = 0: sumr2 = 0: sumd = 0: sumd2 = 0
FOR i% = 1 TO fleetsize%
	correlx = correlx + i%
	correlx2 = correlx2 + i% * i%
NEXT i%
correlx2 = (correlx2 - (correlx ^ 2) / fleetsize%) / (fleetsize% - 1)
IF correlx2 <= 0 THEN
	correlx2 = 0
ELSE
	correlx2 = SQR(correlx2)
END IF
correlx = correlx / fleetsize%

' Ready to run one or more event simulations
eventno% = 1
finished% = 0
DO
PRINT "Event"; eventno%; "races: ";

' Initialise the list of jib numbers
' Boat numbers run from 1 to N, where N is the fleet size
' This is for convenience of evaluating the simulation outcomes
' It is easy to see whether the jib number aligns with the boat's
'  position in the fleet after a simulated race
FOR i% = 1 TO fleetsize%
	fleetjib%(i%) = i%
NEXT i%

' pointslist = Total points scored
' posnsum = Becomes average fleet position over the event
' posnsum2 = Becomes fleet position standard deviation
'   For 80 boats, 5 fleets, 4 promoted, past championship events
'   indicate this should be around 10 or 12 mid-fleet
' resarray = Points scored per race per boat

FOR i% = 1 TO fleetsize%
	pointslist&(i%) = 0
	posnsum#(i%) = 0
	posnsum2#(i%) = 0
	FOR j% = 1 TO maxnumraces%
		resarray%(i%, j%) = 0
	NEXT j%
NEXT i%

' rarray = Correlation coeffs 'r'
' darray = Discrepancy values 'd'
'  where 'c' - calculated over all the races from R1 to Rn
'              ie, given the overall standings in the event after n races
'        'i' - recorded for race Rn alone
'              ie, given the standings in this race only

FOR i% = 1 TO maxnumraces%
	rarrayc(i%) = 0
	rarrayi(i%) = 0
	darrayc(i%) = 0
	darrayi(i%) = 0
NEXT i%

' Assign boats to the seeding races
' You can change seeding assignment by modifying this subroutine
CALL allocseed
FOR i% = 1 TO fleetsize%
	fleetjib%(i%) = newfleet%(i%)
NEXT i%

' Run the seeding races
IF progress% > 0 AND debug% < 1 THEN PRINT "1  ";
CALL raceres1
CALL standings(1)
IF debug% < 0 THEN INPUT a$

IF numraces% > 1 THEN
' Run the second race
IF progress% > 0 AND debug% < 1 THEN PRINT "2 ";
CALL raceres2
CALL standings(2)
IF debug% < 0 THEN INPUT a$
END IF

' Run subsequent races
' There is provision for the event to terminate early if the
' "discrepancy limit" is reached.
IF numraces% > 2 THEN raceno% = 3 ELSE raceno% = numraces% + 1
DO UNTIL raceno% > numraces% OR finished% = 1
	IF progress% > 0 AND debug% < 1 THEN PRINT raceno%;
	CALL raceres3(raceno%)
	CALL standings(raceno%)
	IF debug% < 0 THEN INPUT a$
	IF discrep < discreplimit THEN finished% = 1
	raceno% = raceno% + 1
LOOP

raceno% = raceno% - 1
PRINT
CALL eventstats(eventno%, raceno%)

eventno% = eventno% + 1

LOOP UNTIL eventno% > numevents%

CALL showoverallstats
INPUT "Done"; a$

STOP

handler:
	PRINT "ERROR is actually: ";
	SELECT CASE ERR
CASE 1
PRINT "NEXT without FOR"
CASE 37
PRINT "Argument-count mismatch"
CASE 2
PRINT "Syntax error"
CASE 38
PRINT "Array NOT defined"
CASE 3
PRINT "RETURN without GOSUB"
CASE 40
PRINT "Variable required"
CASE 4
PRINT "Out of DATA"
CASE 50
PRINT "FIELD overflow"
CASE 5
PRINT "Illegal function call"
CASE 51
PRINT "Internal error"
CASE 6
PRINT "Overflow"
CASE 52
PRINT "Bad file name or number"
CASE 7
PRINT "Out of memory"
CASE 53
PRINT "File NOT found"
CASE 8
PRINT "Label NOT defined"
CASE 54
PRINT "Bad file mode"
CASE 9
PRINT "Subscript out of range"
CASE 55
PRINT "File already open"
CASE 10
PRINT "Duplicate definition"
CASE 56
PRINT "FIELD statement active"
CASE 11
PRINT "Division by zero"
CASE 57
PRINT "Device I/O error"
CASE 12
PRINT "Illegal in direct mode"
CASE 58
PRINT "File already exists"
CASE 13
PRINT "TYPE mismatch"
CASE 59
PRINT "Bad record length"
CASE 14
PRINT "Out of string space"
CASE 61
PRINT "Disk full"
CASE 16
PRINT "String formula too complex"
CASE 62
PRINT "Input past end of file"
CASE 17
PRINT "Cannot continue"
CASE 63
PRINT "Bad record number"
CASE 18
PRINT "Function not defined"
CASE 64
PRINT "Bad file name"
CASE 19
PRINT "No RESUME"
CASE 67
PRINT "Too many files"
CASE 20
PRINT "RESUME without error"
CASE 68
PRINT "Device unavailable"
CASE 24
PRINT "Device timeout"
CASE 69
PRINT "Communication-buffer overflow"
CASE 25
PRINT "Device fault"
CASE 70
PRINT "Permission denied"
CASE 26
PRINT "FOR without NEXT"
CASE 71
PRINT "Disk NOT ready"
CASE 27
PRINT "Out of paper"
CASE 72
PRINT "Disk-media error"
CASE 29
PRINT "WHILE without WEND"
CASE 73
PRINT "Feature unavailable"
CASE 30
PRINT "WEND without WHILE"
CASE 74
PRINT "Rename across disks"
CASE 33
PRINT "Duplicate Label"
CASE 75
PRINT "Path/File access error"
CASE 35
PRINT "Subprogram NOT defined"
CASE 76
PRINT "Path NOT found"
END SELECT

INPUT "... <CR> to continue   <CTRL-C> to quit ..."; x$

STOP
RESUME NEXT

END

SUB allocseed
' A simple allocation process is used here, strictly by jib number/
' ability score ranking.

DIM accumposn%(numheats%)

accum% = 1
FOR heatswitch% = 1 TO numheats%
	accumposn%(heatswitch%) = accum%
	accum% = accum% + heatsize1%(heatswitch%)
NEXT heatswitch%

' IF debug% > 0 THEN
'        PRINT "Accum posn: ";
'        FOR i% = 1 TO numheats%
'                PRINT accumposn%(i%);
'        NEXT i%
'        PRINT
' END IF

heatswitch% = 1
posn% = 0
FOR i% = 1 TO fleetsize%
	DO UNTIL posn% < heatsize1%(heatswitch%)
		heatswitch% = heatswitch% + 1
	LOOP
	index% = accumposn%(heatswitch%) + posn%
	newfleet%(index%) = fleetjib%(i%)
	heatswitch% = heatswitch% + 1
	IF heatswitch% > numheats% THEN
		heatswitch% = 1
		posn% = posn% + 1
	END IF
NEXT i%

' IF debug% > 0 THEN
'        PRINT "Seeding allocation: ";
'        FOR i% = 1 TO fleetsize%
'                PRINT newfleet%(i%);
'        NEXT i%
'        PRINT
'END IF

END SUB

SUB calcdiscardn (nr%)
' Calculate number of discards depending upon the number of races
' Note that this calculation needs to be done for each event separately,
'  in case the simulation was terminated 'early' due to low discrepancy
'  score.
' Discard rates from 2 upwards are supported.
' If rate is 2 or 3, then the standard discard at 4 is ignored
'  and discards accumulate after every 2nd or 3rd race.
' If the rate is 4, 5, or 6, the standard discard at 10 is ignored
'  while the standard discard at 4 is maintained.  Discards accumulate
'  every 4th, 5th, or 6th race after race 4.
' If the rate X is 7 or higher, standard discards at 4 and 10 are taken,
'  and then every Xth race after race 10.
' If the "rate" is 0 or 1, then this means there are never any discards (0),
'  or there is only ever one discard (1).
				       
SELECT CASE discrate%
	CASE 0
		numdisc% = 0
	CASE 1
		numdisc% = 1
	CASE IS <= 4
		SELECT CASE nr%
			CASE IS < discrate%
				numdisc% = 0
			CASE IS <= maxnumraces%
				numdisc% = INT(nr% / discrate%)
			CASE ELSE
				numdisc% = INT(maxnumraces% / discrate%)
		END SELECT
	CASE IS <= 6
		SELECT CASE nr%
			CASE IS < 4
				numdisc% = 0
			CASE IS <= maxnumraces%
				numdisc% = 1 + INT((nr% - 4) / discrate%)
			CASE ELSE
				numdisc% = 1 + INT((maxnumraces% - 4) / discrate%)
		END SELECT
	CASE ELSE
		SELECT CASE nr%
			CASE IS < 4
				numdisc% = 0
			CASE IS < 10
				numdisc% = 1
			CASE IS <= maxnumraces%
				numdisc% = 2 + INT((nr% - 10) / discrate%)
			CASE ELSE
				numdisc% = 2 + INT((maxnumraces% - 10) / discrate%)
		END SELECT
END SELECT

IF debug% = -3 THEN PRINT "Number of discards = "; numdisc%

END SUB

SUB calcdiscards (raceno%)

FOR jib% = 1 TO fleetsize%
	IF debug% = -3 THEN PRINT "Discards for "; jib%; ":";
	deduct% = 0
	FOR d% = 1 TO numdisc%
		max% = 0: index% = 0
		FOR r% = 1 TO raceno%
			m% = resarray%(jib%, r%)
			IF m% > max% THEN
				max% = m%
				index% = r%
			END IF
		NEXT r%
		IF debug% = -3 THEN PRINT max%;
		deduct% = deduct% + max%
	      
' Note that each discarded score is set negative in the results array
		resarray%(jib%, index%) = -max%
	NEXT d%
	IF debug% = -3 THEN PRINT " totaling "; deduct%;
      
	IF raceno% <= maxnumraces% THEN
		pointslist&(jib%) = pointslist&(jib%) - deduct%
	ELSE
		IF debug% = -3 THEN PRINT "... but not actually deducted";
	END IF
	IF debug% = -3 THEN INPUT a$
NEXT jib%

END SUB

SUB calcdr (outcome%(), order%, disc, corr)

' Calculate discrepancy and correlation statistics
' order% = 0 if leading boats at the start of the array
'          1 if at the end

' IF debug% > 0 THEN PRINT "Discrepancy values (/place jib disc/): /";
corr = 0
disc = 0
FOR i% = 1 TO fleetsize%
	IF order% = 0 THEN
		posn% = i%
	ELSE
		posn% = fleetsize% - i% + 1
	END IF
	jib% = outcome%(i%)
	corr = corr + jib% * posn%
	d = dwts(posn%) ^ wtpower * ABS(posn% - jib%) ^ dpower
	disc = disc + d
	' IF debug% > 0 THEN PRINT posn%; jib%; d; " / ";
NEXT i%
' IF debug% > 0 THEN INPUT a$

corr = (corr - fleetsize% * correlx * correlx) / (fleetsize% - 1)
corr = corr / (correlx2 * correlx2)

' Normalise the discrepancy value so that it indicates the kind
' of 'average' discrepancy within the fleet between where boats
' have placed in the simulation and where they should be placed
' according to their jib no.
disc = fleetsize% * disc / (fleetsize% ^ wtpower)
disc = disc ^ (1 / dpower)

END SUB

SUB eventstats (eventno%, raceno%)

DIM place(maxfleetsize%), jibno%(maxfleetsize%), outcome%(maxfleetsize%)

' The 'eventno' event is complete after 'raceno' races
' Calculate mean and stdev of placings
' The raw data has been accumulated during the races
IF raceno% > 1 THEN
	FOR i% = 1 TO fleetsize%
	posnsum2#(i%) = (posnsum2#(i%) - ((posnsum#(i%) ^ 2) / raceno%)) / (raceno% - 1)
	IF posnsum2#(i%) <= 0 THEN
		posnsum2#(i%) = 0
	ELSE
		posnsum2#(i%) = SQR(posnsum2#(i%))
	END IF
	posnsum#(i%) = posnsum#(i%) / raceno%
	NEXT i%
ELSE
	FOR i% = 1 TO fleetsize%
		posnsum2#(i%) = 0
	NEXT i%
END IF

' Calculate discards
' Find out how many discards are due, then deduct them
CALL calcdiscardn(raceno%)
IF numdisc% > 0 THEN CALL calcdiscards(raceno%)

IF showlist% = 1 THEN
	IF numdisc% > 0 THEN
		PRINT "Points and positions after"; numdisc%; " discards"
	ELSE
		PRINT "No discards in points and positions"
	END IF
	PRINT "(Av place & StDev include race 1 and ignore any discards)"
	PRINT "Pos   Jib#   Points  Av place   Stdev"
END IF

' Sort the fleet into place order
FOR i% = 1 TO fleetsize%
	place(i%) = pointslist&(i%)
	jibno%(i%) = i%
NEXT i%
CALL sort(fleetsize%, jibno%(), outcome%(), place())

' Calculate discrepancy and agreement statistics
' Leading boats at the end of the array
order% = 1
CALL calcdr(outcome%(), order%, discrep, correlxy)

IF showlist% = 1 THEN
	count% = 1
	FOR i% = fleetsize% TO 1 STEP -1
		PRINT USING "####"; fleetsize% - i% + 1;
		jib% = outcome%(i%)
		PRINT USING "######"; jib%;
		PRINT USING "#########"; pointslist&(jib%);
		PRINT USING "#######.#"; posnsum#(jib%);
		PRINT USING "#######.#"; posnsum2#(jib%)
		count% = count% + 1
		IF count% > 16 THEN
			count% = 1
			INPUT a$
		END IF
	NEXT i%
END IF

PRINT "Final event correlation r =";
PRINT USING "##.###"; correlxy;

' Accumulate overall r & d stats for the simulation
sumr = sumr + correlxy
sumr2 = sumr2 + correlxy ^ 2

PRINT ", final event discrep d =";
PRINT USING "######.#"; discrep
sumd = sumd + discrep
sumd2 = sumd2 + discrep ^ 2

'  Write a CSV output file for the event
IF makefile% = 1 THEN CALL writeeventdatafile(raceno%, outcome%(), eventno%)

IF eventno% <= maxnumevents% THEN
	' Note the event stats per event
	' Pick up the r and d stats without discards from the last race
	' of the event
	' The current r and d stats are calculated after discards
	evarray(eventno%, 1) = rarrayc(raceno%)
	evarray(eventno%, 2) = correlxy
	evarray(eventno%, 3) = darrayc(raceno%)
	evarray(eventno%, 4) = discrep
END IF

IF showlist% = 1 THEN INPUT "Continue"; a$
END SUB

SUB heatres (heatno%, n%, jibno%(), outcome%())
' Simulate a heat

DIM place(maxheatsize%), temp(maxheatsize%)

IF debug% = -2 THEN
	PRINT "Boats in heat"
	FOR i% = 1 TO n%
		PRINT jibno%(i%);
	NEXT i%
	PRINT
END IF

' Get an 'expected' finish score for each boat.
' Higher ability boats get higher scores here.
FOR i% = 1 TO n%
	ability = ablist(jibno%(i%))
	place(i%) = ability + posnsd * (randnorm - .5)
NEXT i%

' IF debug% = -2 THEN
'         PRINT "Expected finishing score (higher = better):"
'         FOR i% = 1 TO n%
'                 PRINT USING "####.#"; place(i%);
'         NEXT i%
'         PRINT
' END IF

' Find expected finishing place for debugging
IF debug% = -2 THEN
	FOR i% = 1 TO n%
		temp(i%) = place(i%)
	NEXT i%
	CALL sort(n%, jibno%(), outcome%(), temp())
	PRINT "Expected finish before incidents: "
	FOR i% = 1 TO n%
		PRINT outcome%(i%);
	NEXT i%
	PRINT
END IF

' Calculate statistical distribution of places
sump = 0    ' Becomes average
sumpp = 0   ' Becomes st.dev.
FOR i% = 1 TO n%
	sump = sump + place(i%)
	sumpp = sumpp + place(i%) ^ 2
NEXT i%
n = n%
sumpp = (sumpp - sump ^ 2 / n) / (n - 1)
sump = sump / n
IF sumpp <= 0 THEN
	sumpp = 0
ELSE
	sumpp = SQR(sumpp)
END IF

' IF debug% > 0 THEN
'         PRINT "Mean="; : PRINT USING "####.#"; sump;
'         PRINT "  StDev="; : PRINT USING "####.#"; sumpp
' END IF

' Generate incidents
numincid% = 0
FOR i% = 1 TO n%
	IF RND < possincid THEN
		numincid% = numincid% + 1
		z = (place(i%) - sump) / sumpp
		incidwt = middlezwt
		IF z < lowerz THEN incidwt = lowerzwt
		IF z > upperz THEN incidwt = upperzwt
		incid = incidwt * (.5 + RND) * incidplaceslost
		IF debug% = -1 OR debug% = -2 THEN
			PRINT "Heat "; heatno%;
			PRINT " INCID to"; jibno%(i%);
			' PRINT ": z=";
			' PRINT USING "###.#"; z;
			' PRINT " with 'raw'=";
			' PRINT USING "###.#"; incid;
			PRINT " losing around";
			pl = n% * incid * incidplaceslost / sumpp
			PRINT USING "###"; pl;
			PRINT " place";
			IF pl = 1 THEN PRINT "  " ELSE PRINT "s  "
		END IF
		' The boat score is decreased by the incident
		place(i%) = place(i%) - incid
	END IF
NEXT i%

' IF debug% = -2 THEN
'         PRINT "Scores after incidents:    "
'         FOR i% = 1 TO n%
'                 PRINT USING "####.#"; place(i%);
'         NEXT i%
'         PRINT
' END IF

' List boats in finishing order
' Boats with high "scores" finish higher
CALL sort(n%, jibno%(), outcome%(), place())

IF debug% = -2 THEN
	PRINT "Finishing order after incidents: "
	FOR i% = 1 TO n%
		PRINT outcome%(i%);
	NEXT i%
	PRINT : INPUT a$
END IF

END SUB

SUB initialparams

' Debug control:
'   0 for normal operation
'  >0 for 'normal' debugging
'   (You will need to 'uncomment' the debug statements you
'   are interested in...)
' Specific values for particular debugs or interests:
'  -1 shows race incidents
'  -2 shows detail of race incident calculations
'  -3 shows discard details
'  -4 shows boat ability list
'  -5 shows race-by-race results
'  -6 shows heat-by-heat results
debug% = 0

' Progress indicator
'  0 for minimal information
'  1 for indicators of where the program has got to
progress% = 1

problem% = 0
DO
	PRINT
	PRINT "Largest permitted heat size (between"; minheatsize%; "and"; maxheatsize%; ") ";
	INPUT largestheatsize%
	IF largestheatsize% < minheatsize% THEN problem% = 1: PRINT "Largest heat size is too small!"
	IF largestheatsize% > maxheatsize% THEN PRINT "Largest heat size is too large!": problem% = 1
	maxnumpromo% = INT(largestheatsize% / 2)
	PRINT "Number of boats to promote/relegate (between"; minnumpromo%; "and"; maxnumpromo%; ") ";
	INPUT numpromo%
	IF numpromo% < 0 THEN PRINT "Don't be silly!": problem% = 1
	IF numpromo% > maxnumpromo% THEN PRINT "Too many promotions!": problem% = 1
	maxfleetsize% = maxnumheats% * (largestheatsize% - numpromo%) + numpromo%
	PRINT "Fleet size (between"; minfleetsize%; "and"; maxfleetsize%; ") ";
	INPUT fleetsize%
	IF fleetsize% < minfleetsize% THEN problem% = 1: PRINT "Fleet too small!"
	IF fleetsize% > maxfleetsize% THEN PRINT "Fleet too large!": problem% = 1
	numheats% = INT((fleetsize% - numpromo%) / (largestheatsize% - numpromo%))
	diff% = fleetsize% - ((numheats% - 1) * (largestheatsize% - numpromo%) + largestheatsize%)
	IF diff% > 0 THEN numheats% = numheats% + 1
	IF numheats% > maxnumheats% THEN PRINT "Too many heats!": problem% = 1
	' IF debug% > 0 THEN PRINT "Num heats ="; numheats%; "  Diff ="; diff%
      
	IF problem% > 0 THEN PRINT "Let's do all this again..."
LOOP UNTIL problem% = 0

PRINT "How many races to simulate in one 'event'";
INPUT numraces%
IF numraces% > maxnumraces% THEN
	PRINT "OK,"; numraces%;
	PRINT "it is, but race-by-race results beyond race"; maxnumraces%;
	PRINT "will not be recorded"
END IF

PRINT "How many 'events' to simulate";
INPUT numevents%
IF numevents% > maxnumevents% THEN
	PRINT "OK,"; numevents%;
	PRINT "it is, but event summary results beyond event"; maxnumevents%;
	PRINT "will not be recorded"
END IF

PRINT "Display final points list after every event (0=No, 1=Yes)";
INPUT showlist%

PRINT "You can create a CSV file for each simulated event, called 'FFFFFxxx.CSV'."
IF numraces% <= maxnumraces% THEN
	PRINT " This contains all"; numraces%; "races of the event"
ELSE
	PRINT " This contains the first"; maxnumraces%; "races of the event"
END IF
PRINT " showing points scored by each boat, and correl r and discrepancy d."
PRINT " The file is 10 KB for every 20 races with 80 boats (it's small!)."
PRINT " Write a race-by-race results file for each event (0=No, 1=Yes)";
INPUT makefile%
IF makefile% = 1 THEN
	PRINT "File name (up to 5 characters) 'FFFFF'";
	INPUT file$
	a$ = file$ + "xxx" + ".CSV"
	PRINT "(Note these files will be called '"; a$; "', where xxx=event number"
	a$ = file$ + "0" + ".CSV"
	PRINT " and the simulation summary file will be called '"; a$; "')"
END IF

IF makefile% = 0 THEN
	PRINT " Write a summary results file for all events (0=No, 1=Yes)";
	INPUT makefilesum%
	IF makefilesum% = 1 THEN
		PRINT "File name (up to 5 characters) 'FFFFF'";
		INPUT file$
		' Eg = "SIM" + LTRIM$(STR$(numpromo%)) + "p"
		a$ = file$ + "0" + ".CSV"
		PRINT "(Note this file will be called '"; a$; "')"
	END IF
END IF

PRINT
PRINT "To change any of the following parameters, stop the program"
PRINT "and reset these values in the 'initialparams' subroutine."
PRINT "------------------"

' Estimate the number of discards.
' Note this subroutine is run again when the event actually ends,
' in case it terminates early and the actual number of discards changes.
discrate% = 9
CALL calcdiscardn(numraces%)
SELECT CASE numdisc%
CASE 0
	disc$ = "There are no discards"
CASE 1
	disc$ = "There is exactly one discard"
CASE ELSE
	disc$ = "There are "
	disc$ = disc$ + LTRIM$(STR$(numdisc%))
	disc$ = disc$ + " discards: after races 4, 10, and every "
	disc$ = disc$ + LTRIM$(STR$(discrate%)) + "th race after"
END SELECT
PRINT disc$

dpower = 2
wtpower = 2
PRINT "Weighting power (suggest between 1.5 and 3) is"; wtpower
PRINT "Discrepancy power (suggest 1 or 2) is"; dpower
PRINT "Discrepancy score, d, is excellently low if <1, good if <";
PRINT USING "##.##"; SQR(fleetsize% ^ (1 / dpower))

discreplimit = 0
IF discreplimit = 0 THEN
	a$ = "The simulation runs all " + LTRIM$(STR$(numraces%))
	a$ = a$ + " races without early termination"
ELSE
	a$ = "The target discrepancy value for early termination is"
	a$ = a$ + LTRIM$(STR$(discreplimit))
END IF
PRINT a$

posnsd = .25
PRINT "Boat ability varies from race to race by"; posnsd * 100; "%"

' Scaling parameter for the random variation in boat ability
abmult = 4
posnsd = abmult * posnsd

' Parameters for race incidents
possincid = .15
lowerz = -.5
lowerzwt = 2
upperz = 1.2
upperzwt = .5
middlezwt = 4
incidplaceslost = .25
PRINT "Incidents (probability="; possincid; "):"
PRINT "  Significance of probable incident when z place in heat > "; upperz; " is"; upperzwt
PRINT "  Significance of probable incident when z place in heat < "; lowerz; " is"; lowerzwt
PRINT "  Significance of probable incident with z between "; lowerz; " and "; upperz; " is"; middlezwt
PRINT "  Average number of places lost in an incident is"; incidplaceslost * 100; "% of heat size"

seed = .4596
PRINT "Random number generator seed (between 0 and 1) is"; seed
RANDOMIZE (seed)

' Set the abilities of the boats
CALL setability
PRINT "Boat ability profile is"; abmethod%; ": '"; abmethod$; "'"
' To change, edit the 'setability' subroutine."

' Set the weights for the discrepancy scores
CALL setdwts
PRINT "Weighting method for boats is"; dwtmethod%; ": '"; dwts$; "'"
' To change, edit the 'setdwts' subroutine."

PRINT "------------------"
INPUT "Ready"; a$
END SUB

SUB makeheatsizes
' heatsize1 is for the seeding race, heatsize2 for the second race,
' and heatsize 3 for the third and subsequent races.

' The following few lines are a quick and dirty way of redoing the
' calculations if it turns out that the Race 2 heat sizes get too big.

numheats% = numheats% - 1
problem% = 0
DO
IF problem% = 1 THEN
	problem% = 0
	PRINT "Hmmm... Increasing number of heats, at least one Race 2 heat is too large..."
END IF

numheats% = numheats% + 1
PRINT "Estimating the number of heats at "; numheats%
PRINT "       ";
FOR i% = 1 TO maxnumheats%
	PRINT "  "; CHR$(64 + i%);
NEXT i%
PRINT

IF numheats% > maxnumheats% THEN
	PRINT "SORRY! Too many heats now!"
	INPUT "Hit <ENTER> to end"; a$
	STOP
END IF

' ------------- Race 1 --------------
' For the seeding race, the fleet is divided equally between the heats,
' with any "odd" boats allocated to the lower heat(s).

num% = INT(fleetsize% / numheats%)
extra% = fleetsize% - num% * numheats%
' IF debug% > 0 THEN PRINT "Seeding:  Num ="; num%; "  Extra ="; extra%
FOR i% = numheats% TO 1 STEP -1
	heatsize1%(i%) = num%
	IF extra% > 0 THEN
		heatsize1%(i%) = heatsize1%(i%) + 1
		extra% = extra% - 1
	END IF
NEXT i%

check% = 0
PRINT "Race 1 ";
FOR i% = 1 TO numheats%
	PRINT USING "###"; heatsize1%(i%);
	check% = check% + heatsize1%(i%)
NEXT i%
IF debug% > 0 THEN PRINT " Check ="; check%;
PRINT

' ------------- Race 2 --------------
' For the second race, the heats are multiples of the number of heats,
' larger multiples in higher heats, and all remaining boats are allocated
' to the lowest heat.

num% = INT((fleetsize% - numpromo%) / numheats%)
num2% = INT(num% / numheats%)
extra% = fleetsize% - num2% * numheats% * numheats%

' IF debug% > 0 THEN PRINT "Second:  Num2 ="; num2%; "  Extra ="; extra%

accum% = 0
FOR i% = 1 TO numheats%
	heatsize2%(i%) = num2% * numheats%
	IF (extra% - numpromo%) >= numheats% THEN
		heatsize2%(i%) = heatsize2%(i%) + numheats%
		extra% = extra% - numheats%
	END IF
	accum% = accum% + heatsize2%(i%)
NEXT i%
heatsize2%(numheats%) = heatsize2%(numheats%) + fleetsize% - accum%

IF heatsize2%(numheats%) > largestheatsize% THEN problem% = 1
IF numheats% > 1 THEN
	FOR i% = 1 TO numheats% - 1
		IF heatsize2%(i%) > (largestheatsize% - numpromo%) THEN problem% = 1
	NEXT i%
END IF

check% = 0
PRINT "Race 2 ";
FOR i% = 1 TO numheats%
	PRINT USING "###"; heatsize2%(i%);
	check% = check% + heatsize2%(i%)
NEXT i%
IF debug% > 0 THEN PRINT " Check ="; check%;
PRINT

' ------------- Race 3+ --------------
' For the third and all following races, the lowest heat is larger by the
' number of promotions, and "odd" boats are allocated starting with the
' higher heat(s).

num% = INT((fleetsize% - numpromo%) / numheats%)
extra% = (fleetsize% - numpromo%) - num% * numheats%
' IF debug% > 0 THEN PRINT "Third plus:  Num ="; num%; "  Extra ="; extra%
FOR i% = 1 TO numheats%
	heatsize3%(i%) = num%
	IF (i% < numheats%) AND (extra% > 0) THEN
		heatsize3%(i%) = heatsize3%(i%) + 1
		extra% = extra% - 1
	END IF
NEXT i%
heatsize3%(numheats%) = heatsize3%(numheats%) + numpromo%

check% = 0
PRINT "Race 3+";
FOR i% = 1 TO numheats%
	PRINT USING "###"; heatsize3%(i%);
	check% = check% + heatsize3%(i%)
NEXT i%
IF debug% > 0 THEN PRINT " Check ="; check%;
PRINT

LOOP UNTIL problem% = 0

END SUB

SUB raceres1
' Simulate the seeding races

DIM jibno%(maxheatsize%), outcome%(maxheatsize%)

FOR heatno% = numheats% TO 1 STEP -1
	heatn% = heatsize1%(heatno%)

	' Find the start and end of the "block" of boats in the heat
	IF heatno% = numheats% THEN
		startplace% = fleetsize% - heatn% + 1
		endplace% = fleetsize%
	ELSE
		startplace% = startplace% - heatn%
		endplace% = startplace% + heatn% - 1
	END IF
       
	IF debug% = -6 THEN
		PRINT
		PRINT "Race 1 Heat "; heatno%; "  n="; heatn%;
		PRINT "  Start="; startplace%; "  End="; endplace%
	END IF

	' Put the boats into a heat list
	FOR i% = startplace% TO endplace%
		jibno%(i% - startplace% + 1) = fleetjib%(i%)
	NEXT i%
       
	IF debug% = -6 THEN
		PRINT "Boats in heat: ";
		FOR i% = 1 TO heatn%
			PRINT jibno%(i%);
		NEXT i%
		PRINT
	END IF

	' Call the heat results simulator
	CALL heatres(heatno%, heatn%, jibno%(), outcome%())
       
	IF debug% = -6 THEN
		PRINT "Heat results:  ";
		FOR i% = 1 TO heatn%
			PRINT outcome%(i%);
		NEXT i%
		PRINT : INPUT a$
	END IF

	' Put the results into the new fleet listing
	' Note that for Race 1 this needs re-arrangement afterwards...
	FOR i% = 1 TO heatn%
		newfleet%(startplace% + i% - 1) = outcome%(i%)
	NEXT i%

NEXT heatno%

' IF debug% > 0 THEN
'         PRINT "Fleet list: ";
'         FOR i% = 1 TO fleetsize%
'                 PRINT newfleet%(i%);
'         NEXT i%
'         PRINT
' END IF

' Allocate the top finishers in each seeding race to the top heat, etc

DIM accumposn%(numheats%)

accum% = 1
FOR heatswitch% = 1 TO numheats%
	accumposn%(heatswitch%) = accum%
	accum% = accum% + heatsize1%(heatswitch%)
NEXT heatswitch%

' IF debug% > 0 THEN
'        PRINT "Accum posn: ";
'        FOR i% = 1 TO numheats%
'                PRINT accumposn%(i%);
'        NEXT i%
'        PRINT
' END IF

' Note that "fleetjib" array will contain the list of jib
' numbers in order of finishing in the race.

heatswitch% = 1
posn% = 0
FOR i% = 1 TO fleetsize%
	DO UNTIL posn% < heatsize1%(heatswitch%)
		heatswitch% = heatswitch% + 1
	LOOP
	index% = accumposn%(heatswitch%) + posn%
	fleetjib%(i%) = newfleet%(index%)
	heatswitch% = heatswitch% + 1
	IF heatswitch% > numheats% THEN
		heatswitch% = 1
		posn% = posn% + 1
	END IF
NEXT i%

IF debug% = -6 THEN
       PRINT "Seeded allocation: ";
       FOR i% = 1 TO fleetsize%
	       PRINT fleetjib%(i%);
       NEXT i%
       PRINT
END IF

' IF debug% > 0 THEN INPUT "Race 1 complete..."; a$

END SUB

SUB raceres2

' Simulate the second race
DIM jibno%(maxheatsize%), outcome%(maxheatsize%)

FOR heatno% = numheats% TO 1 STEP -1
	heatn% = heatsize2%(heatno%)

	' Find the start and end of the "block" of boats in the heat
	IF heatno% = numheats% THEN
		startplace% = fleetsize% - heatn% + 1
		endplace% = fleetsize%
	ELSE
		startplace% = startplace% - heatn%
		endplace% = startplace% + heatn% - 1
	END IF
       
	IF debug% = -6 THEN
		PRINT
		PRINT "Race 2 Heat "; heatno%; "  n="; heatn%;
		PRINT "  Start="; startplace%; "  End="; endplace%
	END IF

	' Put the boats into a heat list
	FOR i% = startplace% TO endplace%
		jibno%(i% - startplace% + 1) = fleetjib%(i%)
	NEXT i%
      

	' Include the promoted boats
	IF heatno% < numheats% THEN
		FOR i% = 1 TO numpromo%
			jibno%(i% + heatn%) = fleetjib%(i% + endplace%)
		NEXT i%
		heatn% = heatn% + numpromo%
	END IF

	IF debug% = -6 THEN
		PRINT "Boats in heat: ";
		FOR i% = 1 TO heatn%
			PRINT jibno%(i%);
		NEXT i%
		PRINT
	END IF


	' Call the heat results simulator
	CALL heatres(heatno%, heatn%, jibno%(), outcome%())
       
	IF debug% = -6 THEN
		PRINT "Heat results:  ";
		FOR i% = 1 TO heatn%
			PRINT outcome%(i%);
		NEXT i%
		PRINT : INPUT a$
	END IF

	' Put the results into the new fleet listing
	' The "fleetjib" array will contain the list of jib
	' numbers in order of finishing in the race.
       
	FOR i% = 1 TO heatn%
		fleetjib%(startplace% + i% - 1) = outcome%(i%)
	NEXT i%

	IF debug% = -6 THEN
		PRINT "Revised fleet: ";
		FOR i% = 1 TO fleetsize%
			PRINT fleetjib%(i%);
		NEXT i%
		PRINT
	END IF

NEXT heatno%

' IF debug% > 0 THEN
'         PRINT "Fleet list: ";
'         FOR i% = 1 TO fleetsize%
'                 PRINT fleetjib%(i%);
'         NEXT i%
'         PRINT
' END IF

' IF debug% > 0 THEN INPUT "Race 2 complete..."; a$

END SUB

SUB raceres3 (racen%)

' Simulate the third and following races
DIM jibno%(maxheatsize%), outcome%(maxheatsize%)

FOR heatno% = numheats% TO 1 STEP -1
	heatn% = heatsize3%(heatno%)

	' Find the start and end of the "block" of boats in the heat
	IF heatno% = numheats% THEN
		startplace% = fleetsize% - heatn% + 1
		endplace% = fleetsize%
	ELSE
		startplace% = startplace% - heatn%
		endplace% = startplace% + heatn% - 1
	END IF
       
	IF debug% = -6 THEN
		PRINT
		PRINT "Race"; racen%; "Heat "; heatno%; "  n="; heatn%;
		PRINT "  Start="; startplace%; "  End="; endplace%
	END IF

	' Put the boats into a heat list
	FOR i% = startplace% TO endplace%
		jibno%(i% - startplace% + 1) = fleetjib%(i%)
	NEXT i%
     

	' Include the promoted boats
	IF heatno% < numheats% THEN
		FOR i% = 1 TO numpromo%
			jibno%(i% + heatn%) = fleetjib%(i% + endplace%)
		NEXT i%
		heatn% = heatn% + numpromo%
	END IF
       
	IF debug% = -6 THEN
		PRINT "Boats in heat: ";
		FOR i% = 1 TO heatn%
			PRINT jibno%(i%);
		NEXT i%
		PRINT
	END IF

	' Call the heat results simulator
	CALL heatres(heatno%, heatn%, jibno%(), outcome%())
       
	IF debug% = -6 THEN
		PRINT "Heat results:  ";
		FOR i% = 1 TO heatn%
			PRINT outcome%(i%);
		NEXT i%
		INPUT a$
	END IF

	' Put the results into the new fleet listing
	' The "fleetjib" array will contain the list of jib
	' numbers in order of finishing in the race.
	FOR i% = 1 TO heatn%
		fleetjib%(startplace% + i% - 1) = outcome%(i%)
	NEXT i%

	IF debug% = -6 THEN
		PRINT "Revised fleet: ";
		FOR i% = 1 TO fleetsize%
			PRINT fleetjib%(i%);
		NEXT i%
		PRINT
	END IF

NEXT heatno%

' IF debug% >0 THEN
'         PRINT "Fleet board: ";
'         FOR i% = 1 TO fleetsize%
'                 PRINT fleetjib%(i%);
'         NEXT i%
'         PRINT : INPUT a$
' END IF

END SUB

FUNCTION randnorm!
' Create a quasi-Gaussian random number
' by simply taking the average of six uniform numbers
x = (RND + RND + RND + RND + RND + RND) / 6
' IF debug% > 0 THEN PRINT x;
randnorm = x
END FUNCTION

SUB setability

' The boat ability should be a number between 1 (highest ability) and
' 0 (lowest ability).

' A number of ability formulae are possible.  Some are given here.
' To chose the one you want to use, simply set the "abmethod" variable.

' 1: closely-spaced linear spread
' 2: widely-spaced linear spread
' 3: hard sharp J
' 4: cosine S
' 5: soft J

abmethod% = 3
SELECT CASE abmethod%
CASE 1
	' Ability setting here is simply the inverse of the jib number
	' Result is a closely-spaced linear spread
	abmethod$ = "Linear, closely spaced"
	FOR i% = 1 TO fleetsize%
		ablist(i%) = (fleetsize% - i%) / maxfleetsize%
	NEXT i%

CASE 2
	' Perhaps small fleets should have relatively more spread abilities.
	' To have the abilities spread, use 'fleetsize%' as the divisor
	' instead of 'maxfleetsize%'.  Result is a more widely-spaced
	' linear spread.
	abmethod$ = "Linear, widely spaced"
	FOR i% = 1 TO fleetsize%
		ablist(i%) = (fleetsize% - i%) / fleetsize%
	NEXT i%

CASE 3
	' An inverse power method yielding a "J" shaped curve
	' Change the early slope of the curve by changing the exponent
	expon = -.25
	abmethod$ = "'Hard' power J, exponent " + LTRIM$(STR$(expon))
	FOR i% = 1 TO fleetsize%
		ablist(i%) = i% ^ expon
	NEXT i%
	lower = ablist(fleetsize%)
	FOR i% = 1 TO fleetsize%
		ablist(i%) = ablist(i%) - lower
	NEXT i%
	upper = ablist(1)
	FOR i% = 1 TO fleetsize%
		ablist(i%) = ablist(i%) / upper
	NEXT i%

CASE 4
	' The first sweep of the Cosine curve
	' Result should be an "S" shaped curve giving two "lumps"
	' of ability, one high and one low.
	abmethod$ = "Cosine S curve"
	FOR i% = 1 TO fleetsize%
		z = 3.141596 * i% / fleetsize%
		ablist(i%) = (1 + COS(z)) / 2
	NEXT i%

CASE 5
	' Another inverse exponential method yielding a "J" shaped curve
	' Change the slope of the curve by changing the exponent
	expon = -.05
	abmethod$ = "'Soft' exponential J, exponent " + LTRIM$(STR$(expon))
	FOR i% = 1 TO fleetsize%
		ablist(i%) = EXP(expon * i%)
	NEXT i%
	lower = ablist(fleetsize%)
	FOR i% = 1 TO fleetsize%
		ablist(i%) = ablist(i%) - lower
	NEXT i%
	upper = ablist(1)
	FOR i% = 1 TO fleetsize%
		ablist(i%) = ablist(i%) / upper
	NEXT i%

END SELECT

IF debug% = -4 THEN
	PRINT "Abilities: ";
	FOR i% = 1 TO fleetsize%
		PRINT USING "#.###"; ablist(i%);
		PRINT "  ";
	NEXT i%
	INPUT a$
END IF

END SUB

SUB setdwts
' Calculate the array of discrepancy weights
' The weights should be between 1 and 'fleetsize'

' 1: Egalitarian (no differential weightings, all the same)
' 2: Very sharp, very steep J  (exp-1)
' 3: Decreasing linear
' 4: Sharp, hard J for top quarter, 0 otherwise (exp-.25)
' 5: Soft J for fleet (exp-.02)

dwtmethod% = 2
SELECT CASE dwtmethod%

CASE 1
	' No weighting for top boats at all.  Everyone the same.
	dwts$ = "Every boat weighted the same"
	FOR i% = 1 TO fleetsize%
		dwts(i%) = fleetsize% / 2
	NEXT i%

CASE 2
	' Exponential weighting for top boats.
	' This is very steep weighting.
	dwts$ = "Top boats very steeply weighted (exp-1)"
	ex = -1!
	FOR i% = 1 TO fleetsize%
		dwts(i%) = i% ^ ex
	NEXT i%
	lower = dwts(fleetsize%)
	FOR i% = 1 TO fleetsize%
		dwts(i%) = dwts(i%) - lower
	NEXT i%
	upper = dwts(1)
	FOR i% = 1 TO fleetsize%
		dwts(i%) = dwts(i%) / upper
		dwts(i%) = fleetsize% * dwts(i%)
	NEXT i%

CASE 3
	' Linear weighting for top boats.
	dwts$ = "Linear weighting for top boats"
	FOR i% = 1 TO fleetsize%
		dwts(i%) = fleetsize% - i% + 1
	NEXT i%

CASE 4
	' Exponential weighting for top boats, other boats zero.
	' Relatively steep weighting.
	ex = -.25
	FOR i% = 1 TO fleetsize%
		dwts(i%) = i% ^ ex
	NEXT i%
	lower = dwts(fleetsize%)
	FOR i% = 1 TO fleetsize%
		dwts(i%) = dwts(i%) - lower
	NEXT i%
	upper = dwts(1)
	FOR i% = 1 TO fleetsize%
		dwts(i%) = dwts(i%) / upper
		dwts(i%) = fleetsize% * dwts(i%)
	NEXT i%
       
	' Now ignore the bottom three-quarters of the fleet
	ig% = fleetsize% / 4
	dwts$ = "Bottom " + LTRIM$(STR$(fleetsize% - ig%)) + " boats ignored,"
	dwts$ = dwts$ + " only top " + LTRIM$(STR$(ig%)) + " boats steeply weighted"
	FOR i% = ig% + 1 TO fleetsize%
		dwts(i%) = 0
	NEXT i%

CASE 5
	' Exponential weighting for top boats.
	' Relatively soft J weighting.
	dwts$ = "Top boats weighted with 'soft' J (exp-.02)"
	ex = -.02
	FOR i% = 1 TO fleetsize%
		dwts(i%) = i% ^ ex
	NEXT i%
	lower = dwts(fleetsize%)
	FOR i% = 1 TO fleetsize%
		dwts(i%) = dwts(i%) - lower
	NEXT i%
	upper = dwts(1)
	FOR i% = 1 TO fleetsize%
		dwts(i%) = dwts(i%) / upper
		dwts(i%) = fleetsize% * dwts(i%)
	NEXT i%
      
END SELECT

END SUB

SUB showoverallstats

' The simulation run is complete
' Calculate the overall simulation statistics

IF numevents% > 1 THEN
sumr2 = (sumr2 - ((sumr ^ 2) / numevents%)) / (numevents% - 1)
IF sumr2 <= 0 THEN
	sumr2 = 0
ELSE
	sumr2 = SQR(sumr2)
END IF
sumd2 = (sumd2 - ((sumd ^ 2) / numevents%)) / (numevents% - 1)
IF sumd2 <= 0 THEN
	sumd2 = 0
ELSE
	sumd2 = SQR(sumd2)
END IF
ELSE
sumr2 = 0: sumd2 = 0
END IF

sumr = sumr / numevents%
sumd = sumd / numevents%

PRINT "After"; numevents%; "events the overall simulation stats are:"
PRINT "Correlation avg = "; sumr; ", st.dev = "; sumr2
PRINT "Discrepancy avg = "; sumd; ", st.dev = "; sumd2

'  Write a CSV output file for the event summary
IF makefile% = 1 OR makefilesum% = 1 THEN CALL writesumdatafile

END SUB

SUB sort (n%, jibno%(), outcome%(), place())
' Take the arbitrary results in the "place" array
' and sort them into order.
' On entry to this subroutine, the jib no of each of the
' boats listed in the "places" array is recorded in the
' "jibno" array.  The jib no of the boat whose score is
' places(2) is jibno(2) for example.
' The result of the sort is to insert the jib nos into
' the outcome array in order of their 'places' score,
' boats with high scores listed first.
' This order is the "natural" order, from high scores to low.
' It is NOT the "sailing" order, where the low scores are best.
' NB The "places" array is destroyed by the sort.
' NB The "jibno" array has no meaning after the sort.

nextpos% = 1
DO
	found% = 0
	max = -10000
	FOR i% = 1 TO n%

		' Find the (next) highest place
		IF place(i%) > max THEN
			found% = 1
			max = place(i%)
			maxi% = i%
		END IF
	NEXT i%
	IF found% > 0 THEN
		outcome%(nextpos%) = jibno%(maxi%)
		' IF debug% > 0 THEN PRINT nextpos%; "-"; outcome%(nextpos%); " ";
		nextpos% = nextpos% + 1
		IF nextpos% > n% + 1 THEN
			PRINT "Houston, we have a problem in the SORT routine..."
			STOP
		END IF
		place(maxi%) = -10001
	END IF
LOOP UNTIL found% = 0

END SUB

SUB standings (raceno%)

DIM place(maxfleetsize%), jibno%(maxfleetsize%), outcome%(maxfleetsize%)

' ----------- Step 1 -------------
' Assign each boat its overall points for the race.
' The race results are held in the 'fleetjib' array
' which is the list of jib numbers in order of finishing
' in the race.
' Accumulate the boat's points, accumulate the raw data for places
' Enter the boat's points for the race into the results 'resarray' array

IF raceno% = 1 THEN
	i% = 1
	p% = 1
	DO
		FOR heatno% = 1 TO numheats%
			IF i% <= fleetsize% THEN
				jib% = fleetjib%(i%)
				pointslist&(jib%) = p%
				resarray%(jib%, 1) = p%
				posnsum#(jib%) = i%
				posnsum2#(jib%) = i% ^ 2
				i% = i% + 1
			END IF
		NEXT heatno%
		p% = p% + 1
	LOOP UNTIL i% > fleetsize%
ELSE
	FOR i% = 1 TO fleetsize%
		jib% = fleetjib%(i%)
		pointslist&(jib%) = pointslist&(jib%) + i%
		IF raceno% <= maxnumraces% THEN
			resarray%(jib%, raceno%) = i%
		END IF
		posnsum#(jib%) = posnsum#(jib%) + i%
		posnsum2#(jib%) = posnsum2#(jib%) + i% ^ 2
	NEXT i%
END IF

' --------------- Step 2 ------------------
' Calculate correlation & discrepancy between jib no and place
' just for this race.  The list of the jib nos is in "fleetjib" array.
' Leading boats at the start of the array
order% = 0
CALL calcdr(fleetjib%(), order%, discrep, correlxy)

IF raceno% <= maxnumraces% THEN
	' Make a note of the race's r and d statistics
	rarrayi(raceno%) = correlxy
	darrayi(raceno%) = discrep
END IF

IF debug% = -5 THEN
	PRINT "For this race, r = "; correlxy;
	PRINT "   discrep = "; discrep;
END IF

' --------------- Step 3 ------------------
' Find out the overall standings of the boats after raceno races.
' Put the total points into the 'place' array, and note the
' jib no corresponding to these total points in the 'jibno' array.
FOR i% = 1 TO fleetsize%
	place(i%) = pointslist&(i%)
	jibno%(i%) = i%
NEXT i%

CALL sort(fleetsize%, jibno%(), outcome%(), place())

' Note that the results come back in "natural" order of points
' scored, from high to low.  This order must be listed in reverse
' for sailing results.
IF debug% = -5 THEN
	PRINT "Results (no discards yet) after race "; raceno%
	FOR i% = fleetsize% TO 1 STEP -1
		PRINT outcome%(i%);
	NEXT i%
	PRINT
END IF

' Calculate correlation & discrepancy between jib no and standing
' for the overall event up to and including this race (no discards!)
' Leading boats at the end of the array
order% = 1
CALL calcdr(outcome%(), order%, discrep, correlxy)

IF raceno% <= maxnumraces% THEN
	' Record the r and d statistics for the event up to this race
	rarrayc(raceno%) = correlxy
	darrayc(raceno%) = discrep
END IF

IF debug% = -5 THEN
	PRINT "Event so far, r = "; correlxy;
	PRINT "   discrep = "; discrep;
	INPUT a$
END IF

END SUB

SUB writeeventdatafile (raceno%, outcome%(), eventno%)

IF raceno% > maxnumraces% THEN raceno% = maxnumraces%

extn$ = ".CSV"
seq$ = LTRIM$(STR$(eventno%))
fileout$ = file$ + seq$ + extn$

OPEN fileout$ FOR OUTPUT AS #2

PRINT #2, "Position, Jib,";
FOR race% = 1 TO raceno%
	PRINT #2, USING "####"; race%;
	PRINT #2, ",";
NEXT race%
PRINT #2, " Total, Avg, StDev"

FOR i% = fleetsize% TO 1 STEP -1
	PRINT #2, USING "###"; fleetsize% - i% + 1;
	PRINT #2, ",";
	jib% = outcome%(i%)
	PRINT #2, USING "###"; jib%;
	PRINT #2, ",";
      
	FOR race% = 1 TO raceno%
		PRINT #2, USING "####"; resarray%(jib%, race%);
		PRINT #2, ",";
	NEXT race%

	PRINT #2, USING "#########"; pointslist&(jib%);
	PRINT #2, ",";
	PRINT #2, USING "#######.#"; posnsum#(jib%);
	PRINT #2, ",";
	PRINT #2, USING "#######.#"; posnsum2#(jib%)
NEXT i%

PRINT #2, "Ind. Correl r, ,";
FOR race% = 1 TO raceno%
	PRINT #2, USING "##.###"; rarrayi(race%);
	IF race% < raceno% THEN PRINT #2, ",";
NEXT race%
PRINT #2, ""
PRINT #2, "Ind. Discrep d, ,";
FOR race% = 1 TO raceno%
	PRINT #2, USING "######.#"; darrayi(race%);
	IF race% < raceno% THEN PRINT #2, ",";
NEXT race%
PRINT #2, " "

PRINT #2, "Cum. Correl r, ,";
FOR race% = 1 TO raceno%
	PRINT #2, USING "##.###"; rarrayc(race%);
	PRINT #2, ",";
NEXT race%

' NOTICE
' While the per-race correlations & discrepancies are calculated
' ignoring any discards, this figure is calculated after discards
' have been taken.
PRINT #2, USING "##.####"; correlxy

PRINT #2, "Cum. Discrep d, ,";
FOR race% = 1 TO raceno%
	PRINT #2, USING "######.#"; darrayc(race%);
	PRINT #2, ",";
NEXT race%

' NOTICE
' While the per-race correlations & discrepancies are calculated
' ignoring any discards, this figure is calculated after discards
' have been taken.
' It is therefore useful to compare it to see the effect of discards.
PRINT #2, USING "######.#"; discrep

CLOSE #2

IF numraces% > maxnumraces% THEN
	PRINT "Note that the event involved"; numraces%; "races,"
	PRINT "but only the first"; raceno%; "race results have been filed."
END IF
PRINT "Data file '"; fileout$; "' written.  ";

END SUB

SUB writesumdatafile

'  Write a CSV output file for the event summary

extn$ = ".CSV"
seq$ = LTRIM$(STR$(0))
fileout$ = file$ + seq$ + extn$

OPEN fileout$ FOR OUTPUT AS #2

PRINT #2, "Event, Correl (no disc), Correl (discs), Discrep (no disc), Discrep (disc)"

ne% = numevents%
IF numevents% > maxnumevents% THEN ne% = maxnumevents%
FOR i% = 1 TO ne%
	PRINT #2, USING "####"; i%;
	PRINT #2, ",";
	PRINT #2, USING "##.###"; evarray(i%, 1);
	PRINT #2, ",";
	PRINT #2, USING "##.###"; evarray(i%, 2);
	PRINT #2, ",";
	PRINT #2, USING "######.#"; evarray(i%, 3);
	PRINT #2, ",";
	PRINT #2, USING "######.#"; evarray(i%, 4)
NEXT i%

IF numevents% > maxnumevents% THEN
	PRINT "Note that"; numevents%; "were simulated,"
	PRINT "but only the first"; ne%; "event summaries have been filed."
END IF

CLOSE #2
PRINT "Data file '"; fileout$; "' written.  ";

END SUB