C THIS IS A SIMULATION WHICH SIMULATES C RATS BAR PRESSING FOR A REINFORCEMENT C OF ELECTRICAL CURRENT APPLIED DIRECTLY C TO THEIR BRAINS VIA CHRONICALLY IMPLANTED C ELECTRODES. C C THIS SUBROUTINE WAS WRITTEN BY MARK KLEIN C AND PROGRAMMED BY MARK KLEIN WITH THE AID C OF RICH NUSSLOCH IN JANUARY OF 1974 C C C THIS IS THE FIRST VERSION SO PLEASE PARDON C THE MISTAKES. SUBROUTINE MODEL(N,NCONDS,*) C C C@@@@@@ SECTION 1 HAS ENTIRELY TO DO WITH HOUSE C KEEPING CHORES AND THE LIKE. C C C VARIABLE TYPE DECLARATION: C INTEGER*4 PHASE,POLES,DIRECT,SHAPE,ELECTR,LCYCLE INTEGER*4 AREA,LESION,DL1,DRUG,EXPER,N,NDED,LIVE INTEGER*4 OPDED1,OPDED2,NBAD,IPROB1,IPROB2,IPROB3 INTEGER*4 MICROC,ITEMP,FIRSTT(11),CLOCK,TTIME INTEGER*4 EXPER1,EXPER2 INTEGER*4 NOV,SINK,IDEV1,IDEV2,IDEV3,IV(12) REAL*4 MICROA,FREQ,PULSE,TRAIN,TESTH,RTEMP,ANDOSE REAL*4 DDOSE,INJECT,FDEPRV,PDEAD,OPER,RATE,T1RATE REAL*4 T2RATE,FIRSTE(11),DARK(24),DIV,OVERD(5) REAL*4 TEFF,DREFF(4,35),SHOT(4),MORT,LMFB(7,7) REAL*4 DMISS(6),DV(6),PRESS(20,2,6),PULL,SEX REAL*4 SIGN(5),FV(12) LOGICAL*4 FLGS(6),EOVER C C COMMON STATEMENTS C COMMON /VARS1/IV,FV/IO1/NOV,SINK/IO/IDEV1,IDEV2,IDEV3 C C DATA STATEMENTS C !! DATA DMISS/6*350./,FLGS/6*.FALSE./ !! DATA EOVER/.FALSE./ C C EQUIVALENCE STATEMENTS C EQUIVALENCE (MICROA,FV(1)),(PHASE,IV(1)),(POLES,IV(2)), 1(DIRECT,IV(3)),(SHAPE,IV(4)),(FREQ,FV(2)),(PULSE,FV(3)), 2(TRAIN,FV(4)),(ELECTR,IV(5)),(LCYCLE,IV(6)),(TESTH,FV(5)), 3(RTEMP,FV(6)),(AREA,IV(7)),(ANDOSE,FV(8)),(LESION,IV(8)), 4(DL1,IV(9)),(DRUG,IV(10)),(DDOSE,FV(9)),(INJECT,FV(10)), 5(FDEPRV,FV(11)),(EXPER1,IV(11)),(SEX,FV(7)) C C ************************MTS SPECIFIC SECTION ******************** C @@@@@@@@@SECTION TO CONTACT CLOCK FOR TIME OF DAY DEFAULT C REMOVE THIS WHOLE SECTION OR MODIFY APPROPRIATELY TO FIT C THE RESIDENT SYSTEM. C INTEGER*4 IT INTEGER*4 JM EQUIVALENCE (IT,JM) C C ******************************END MTS SECTION ******************* C C SET UP FOR DATA FILES C DOUBLE PRECISION MODNM1,MODNM2,MODNM3 COMMON /DATANM/ MODNM1,MODNM2,MODNM3 DATA MODNM1/'RATS.DAT'/,MODNM2/'RATS3.DAT'/,MODNM3/'RATS.BIN'/ DATA DMISS/6*350./,FLGS/6*.FALSE./ DATA EOVER/.FALSE./ C NCONDS=NCONDS! A DEC-10 BUG MAKES THIS NECESSARY C C NOW A LITTLE ERROR CHECKING!!! C 1 IF (SHAPE.EQ.1) GO TO 5 IF (PULSE/1000.LT.1/FREQ) GO TO 5 WRITE(SINK,800) RETURN C@@@@@@@@SECTION 2 C THIS SECTION MAKES THE STUDENT AWARE C OF SOME OF THE TYPICAL PROBLEMS THAT C CAN ARRISE IN THE OPERATION DURING C WHICH THE ELECTRODES ARE IMPLANTED. C C 5 NDED=0 C SECTION 2.1 C THE ANESTHETIC DOSE???? IF (ANDOSE.LE.50) GO TO 10 PDEAD=(ANDOSE-50)*0.05 IF (PDEAD.GT.1.) PDEAD=1. C C FIND THE NUMBER OF ANIMALS TO DIE FRON O.D. C CALL BINOM(N,PDEAD,NDED) NDED=MIN0(NDED,N) C LET THE PEOPLE KNOW HOW MANY DIED. IF (NDED.GT.0) WRITE(SINK,801)NDED C C IS THE DOSE TOO STRONG TO CONTINUE C IF (NDED.LT.N/2) GO TO 11 WRITE(SINK,802) RETURN C ANDOSE TOO LOW? 10 IF(ANDOSE.GT.30)GO TO 11 WRITE(SINK,816) RETURN 11 LIVE=N-NDED C C EXPER=EXPER1+EXPER2 OPER=1-EXPER/1000 OPER=AMAX1(.015,OPER) OPER=AMIN1(.55,OPER) C THIS PART OF THE SECTION INCORPORATES C EXPERIENCE AND CHANCE TO FIND OUT HOW C MANY OF THE LIVING ANIMALS WILL MAKE C IT THROUGH THE OPERATION. CALL BINOM(LIVE,OPER,NDED) LIVE=LIVE-NDED C NOW DETERMINE HOW MANY SUFFER FROM THE C FIRST MALADY CALL URAND2(0,NDED,OPDED1) NDED=NDED-OPDED1 C C INFORM THE FOLKS C IF(OPDED1.GT.0)WRITE(SINK,803)OPDED1 IF(NDED.EQ.0)GO TO 50 C THE NEXT MALADY!!!!!! CALL URAND2(0,NDED,OPDED2) NDED=NDED-OPDED2 C INFORM THE FOLKS!!!!!! IF(OPDED2.GT.0)WRITE(SINK,804)OPDED2 C THE LAST MALADY!!! IF (NDED.GT.0) WRITE(SINK,805)NDED C C C SECTION 2.3 C THIS SECTION SIMULATES SOME OF THE VARIOUS C PROBLEMS THAT THE EXPERIMENTER MIGHT ENCOUNTER C ONCE HE IS ABLE TO BEGIN TRAINING THE RATS C AND AFTERWARDS WHEN HE IS VARIFYING THE C LOCATION OF THE ELECTRODES 50 OPER=OPER/5 C HOW MANY RATS WILL BE EFFECTED?????? CALL BINOM(LIVE,OPER,NBAD) IF(NBAD.EQ.0) GO TO 100 LIVE=LIVE-NBAD C HOW MANY RATS HAVE THE FIRST PROBLEM??? CALL URAND2(0,NBAD,IPROB1) NBAD=NBAD-IPROB1 C LET THE FOLKS KNOW. IF(IPROB1.GT.0) WRITE(SINK,806)IPROB1 IF(NBAD.EQ.0) GO TO 100 C THE SECOND MALADY CALL URAND2(0,NBAD,IPROB2) NBAD=NBAD-IPROB2 C LET THE FOLKS KNOW IF(IPROB2.GT.0) WRITE(SINK,807)IPROB2 IF(NBAD.EQ.0) GO TO 100 C THE THIRD MALADY CALL URAND2(0,NBAD,IPROB3) NBAD=NBAD-IPROB3 C LET THE FOLKS KNOW IF (IPROB3.GT.0)WRITE(SINK,808)IPROB3 C THE LAST ONE IF(NBAD.GT.0)WRITE(SINK,809)NBAD C LET THE FOLKS KNOW HOW MANY THEY HAVE LEFT TO C WORK WITH 100 WRITE(SINK,810)LIVE C C IN THE EVENT THAT ALL OF THE SUBJECTS WERE C WIPED OUT THERE IS NO SENSE IN GOING ANY FURTHER IF(LIVE.EQ.0)RETURN C C C@@@@@@SECTION 3A C THIS SECTION PROCESSES THE FIRST PART OF THE WORK C RELATING TO PHASE,POLES,AND DIRECTION. IF(PHASE.EQ.2) WRITE(SINK,811) IF(POLES.EQ.2) GO TO 200 IF((PHASE.EQ.1).AND.(DIRECT.EQ.2)) MICROA=MICROA-150 IF(MICROA.LT.0)MICROA=1. C C C@@@@@@SECTION 4 C THIS SECTION CALCULATES THE GROSS ENERGY INPUT USED AS C STIMULUS TO THE RAT EACH TIME IT PRESSES THE BAR. C BASED ON THEIR MICROCOULOMB VALUE, THE BASE RATE WILL BE C OBTAINED FROM THE BASE RATE DATA MATRIX. 200 IF(TRAIN.GE.5.) GO TO 210 GO TO 220 C STIMULUS IS SO LONG THAT IT BECOMES AVERSIVE 210 RATE=.001 GO TO 300 C MICROCOULOMBS BY RECTANGULAR CURRENT 220 MICROC=MICROA*FREQ*TRAIN*PULSE*0.001 221 IF(MICROC.LT.21) GO TO 222 GO TO 210 C GET BASE RATE FROM THE MATRIX 222 ITEMP=MICROC IF(AREA.EQ.7) AREA=6 IF(ITEMP.LE.0) GO TO 210 RATE= PRESS(ITEMP,1,AREA) C AN EFFECT DUE JUST TO USING SINE WAVES 255 IF(SHAPE.EQ.1) RATE=RATE*1.5 C C C@@@@@@SECTION 3B C FURTHER CALCULATION OF THE EFFECT OF POLES AND PHASE IF(POLES.EQ.1.AND.PHASE.EQ.2) RATE=RATE*.9 C C C C@@@@@@SECTION 5 C THIS SECTION CALCULATES THE TIME OF DAY EFFECT BASED ON C THE SPECIFIED LIGHT CYCLE, THE TEST HOUR, AND THE INTENSITY C OF THE STIMULATION. C C C INITIALIZE TEMPORARY RATE FOR USE IN FIGURING THE EFFECT OF C THE LIGHT CYCLE 300 T1RATE= -109 ITEMP=(MICROA/100.)+1. C SECTION 5.1 C CHECK TO SEE IF RATE COMES UNDER THE EFFECT DUE TO INITIAL C EXPERIENCE IN THE APPARATUS OR AFTERWARDS WHEN THE RATE HAS C SETTLED DOWN TO A CYCLIC SORT OF EVENT IF(FIRSTT(ITEMP).LT.TESTH) GO TO 315 C IF TESTING OCURRS DURING THAT INITIAL DECLINE T1RATE=RATE*(1.25-FIRSTE(ITEMP)*TESTH**2/100) C C IS THE LIGHT CYCLE SET FOR A CONSTANT LIGHT LEVEL? C 315 IF(LCYCLE.EQ.13) GO TO 350 IF(LCYCLE.EQ.14) GO TO 350 C A TWELVE HOUR LIGHT CYCLE SPECIED? IF (LCYCLE.LT.13) GO TO 325 C LIGHT CYCLE DEFAULT C ********************** MTS SPECIFIC SECTION ******************** C SECTION TO FIND THE TIME OF DAY IN HOURS--USES MTS SPECIFIC C SUBROUTINE "TIME" TO FILL IN THE IT ARRAY WITH THE TIME OF C DAY. C REMOVE THIS SECTION OR MODIFY APPROPRIATELY FOR THE RESIDENT C SYSTEM. C CALL TIME(HRS) DECODE(5,320,HRS) IT 320 FORMAT(I2,3X) C**************END MTS SPECIFIC SECTION********************* CLOCK=JM-1 IF(CLOCK.GT.6) LCYCLE=-CLOCK+19 IF (CLOCK.LE.6) LCYCLE=-CLOCK-5 C RESCALE TESTH TO CONFORM TO A STANDARD 24 HOUR SCHEME 325 TTIME=TESTH+12-LCYCLE TTIME=TTIME-(TTIME/24)*24 C IS THE TEST HOUR DURING THE DARK HALF OF THE CYCLE? IF(TTIME.GE.12) GO TO 375 GO TO 360 C RESCALE TESTH FOR CONSTANT LIGHT LEVELS 350 TTIME= TESTH-(TESTH/24)*24 IF(LCYCLE.EQ.14) GO TO 375 C RESCALE MICROA TO DETERMINE THE INTERACTION 360 ITEMP=MICROA/250.+1. TTIME=(TTIME*15)+180 C CALCULATE THE CIRCADIAN EFFECT T2RATE=RATE*(1+SIGN(ITEMP)*SIN(TTIME/57.29576)) GO TO 390 C FIND THE EFFECT FOR THE HOUR OF DARKNESS 375 T2RATE= RATE*DARK(TTIME+1) C DECIDE WHICH EFFECT TO USE: CHOOSE THE LARGER OF THE TWO 390 RATE=T1RATE IF(T1RATE.LT.T2RATE) RATE=T2RATE C C C C@@@@@@SECTION 6 C THIS SECTION CALCULATES THE EFFECT DUE TO ELECTRODE MATERIAL C DIFFERENCES IF(ELECTR.EQ.2) GO TO 395 IF (TESTH.LE.19) GO TO 395 C INTERACTION OF PHASE AND ELECTRODE MATERIAL IF(PHASE.EQ.2) RATE=RATE*.80 C INTERACTION OF PHASE, DIRECTION, AND MATERIAL IF(PHASE.EQ.1.AND.DIRECT.EQ.2) RATE=RATE*.80 C INTERACTION OF PHASE, POLES, AND MATERIAL IF((POLES.EQ.1).AND.(PHASE.EQ.2)) RATE=RATE*.90 C C C@@@@@@ SECTION 7 C THIS SECTION HANDLES THE EFFECT DUE TO THE AMBIENT ROOM TEMPE 395 CONTINUE 396 IF((RTEMP.GT.22.).AND.(RTEMP.LT.28)) GO TO 400 C IS THE TEMPERATURE TOO LOW? IF(RTEMP.GE.2) GO TO 398 WRITE(SINK,813) RETURN 398 IF(AREA.GT.5) GO TO 400 C RESCALE TEMPERATURE JUST IN CASE THERE WASN'T ENOUGH TIME C FOR THE RAT TO COOL OFF OR WARM UP COMPLETELY IF(TESTH.LT.1) RTEMP=RTEMP+(25-RTEMP)/2 IF(RTEMP.LT.25) GO TO 399 C CALCULATE THE EFFECT DUE TO A HIGH ROOM TEMP RATE=RATE*(1-(RTEMP**2-30*RTEMP+225)/300.0) GO TO 400 C CALCULATE THE EFFECT DUE TO A LOW ROOM TEMP 399 DIV=833.33 IF((AREA.EQ.2).OR.(AREA.EQ.3)) DIV=2500.0 RATE=RATE*(1-(RTEMP**2-50*RTEMP+625)/DIV) C C C C@@@@@@SECTION 8 C THIS SECTION DEALS WITH THE DRUG VARIABLE 400 IF(DDOSE.EQ.0) GO TO 500 IF(DRUG.EQ.6) GO TO 500 C ARE THERE ANY OVERDOSES OF THE PARTICULAR DRUG? IF(DDOSE.LE.OVERD(DRUG)) GO TO 410 C YES??????? PDEAD=(DDOSE-OVERD(DRUG)/(OVERD(DRUG)/2)) IF(PDEAD.GT.1) PDEAD=1. C DETERMINE THE NUMBER OF SUBJECTS TO DIE BY THAT OVERDOSE CALL BINOM(LIVE,PDEAD,NDED) IF (NDED.EQ.0) GO TO 410 LIVE=LIVE-NDED WRITE(SINK,814)NDED IF(LIVE.EQ.0) RETURN C SALINE? 410 IF(DRUG.EQ.5) GO TO 500 C FIGURE TIME DELAY UNTIL TESTING TIM= TESTH-INJECT C THAT TIME DELAY TOO LONG THE DRUGS NO LONGER EFFECTIVE???? IF(TIM.GT.15) GO TO 500 IF(TIM.GT.2) GO TO 415 C EFFECT OF TIME ON THE ACTION OF THE DRUG? TEFF=TIM/2 GO TO 416 415 TEFF=1-(TIM-2)/13 C RESCALE DRUG DOSE SO AS TO USE A MATRIX 416 IDOSE=DDOSE*SHOT(DRUG) 460 RATE=RATE*(1+TEFF*DREFF(DRUG,IDOSE)) C C C C@@@@@@SECTION 9 C THIS SECTION PROCESSES THE LESION VARIABLE 500 IF(LESION.EQ.8) GO TO 600 C SET UP THE MORTALITY RATE FOR THE LESION OPERATION MORT=.20 IF(DL1.NE.08) MORT=.5 CALL BINOM(LIVE,MORT,NDED) IF (NDED.GT.0) WRITE(SINK,815)NDED LIVE=LIVE-NDED C DID THEY WAIT TOO LONG TO TEST FOR THE EFFECT OF THE LESION?? IF(TESTH.LT.168) GO TO 520 IF((LESION.EQ.1.AND.AREA.EQ.6).OR.(DL1.EQ.1.AND.AREA.EQ.6))GOTO510 IF((LESION.EQ.2.AND.AREA.EQ.1).OR.(DL1.EQ.2.AND.AREA.EQ.1))GOTO510 IF((LESION.EQ.4.AND.AREA.EQ.3).OR.(DL1.EQ.4.AND.AREA.EQ.3))GOTO510 GO TO 600 510 RATE=RATE*.5 GO TO 600 520 IF((LESION.EQ.1.AND.AREA.EQ.6).OR.(DL1.EQ.1.AND.AREA.EQ.6))GOTO590 IF((LESION.EQ.2.AND.AREA.EQ.1).OR.(DL1.EQ.2.AND.AREA.EQ.1))GOTO590 IF((LESION.EQ.4.AND.AREA.EQ.3).OR.(DL1.EQ.4.AND.AREA.EQ.3))GOTO590 GO TO(530,540,600,530,600,530,600),AREA C MFB AND RELATED AREA LESIONS 530 RATE=RATE*LMFB(LESION,DL1) GO TO 600 C SEPTAL AREA LESIONS 540 IF((LESION.EQ.4) .OR.(DL1.EQ.4)) RATE=RATE*.5 GO TO 600 C IF THE TEST AREA IS THE SAME ONE THAT THEY DID THE LESION IN. 590 RATE=RATE*.1 C C CCCCCCCSECTION 10 C THIS SECTION HANDLES THE FOOD DEPREVATION VARIABLE 600 IF(FDEPRV.LE.1) GO TO 700 C CALCULATE THE EFFECT IF(FDEPRV.GE.36)RATE=RATE*(.25*FDEPRV-4.4) IF((FDEPRV.GE.20.AND.FDEPRV.LT.36)) FDEPRV=20. IF(FDEPRV.LE.20.)RATE=RATE*SQRT(FDEPRV) 700 IF(RATE.GT.350) RATE=355 701 ITEMP=MICROC IF(ITEMP.EQ.0) ITEMP=1 FLGS(1)=.TRUE. CALL FPOUT1(FLGS,DMISS,&999) C SIMULATE THE DATA FOR THE RATS THAT ACTUALLY MADE IT C DO 799 K=1,LIVE CALL NRAND(RATE,PRESS(ITEMP,2,AREA),PULL) IF(PULL.LT.0.05) PULL=ABS(PULL) IF(PULL.GT.350) PULL=350 EOVER=(PULL.EQ.350).OR.EOVER DV(1)=PULL C OUTPUT THAT SCORE!! KTMP=K CALL FPOUT2(KTMP,1,FLGS,DV,&999) 799 CONTINUE CALL FPOUT3(1) IF(EOVER) WRITE(SINK,817) C INCREMENT THE EXPERIENCE FACTOR C 7999 EXPER2=EXPER2+N EOVER=.FALSE. 800 FORMAT(' IMPOSSIBLE PULSE AND FREQUENCY COMBINATION') 801 FORMAT(' ',I3,' ANIMALS DIED FROM AN ANESTHETIC OVERDOSE') 802 FORMAT(' KILLER ANESTHETIC DOSE TERMINATE EXPERIMENT') 803 FORMAT(' ',I3,' ANIMALS DIED:YOU SEVERED THE MID-SAGGITAL' * ' ARTERY.') 804 FORMAT(' ',I3,' ANIMALS DIED: YOU LET BODY TEMP GET TOO LOW.') 805 FORMAT(' ',I3,' ANIMALS JUST DIED: REASON UNKNOWN.') 806 FORMAT(' ',I3,' ANIMALS NO DATA: ELECTRODES KNOCKED OUT.') 807 FORMAT(' ',I3,' ANIMALS NO DATA: MISSED BRAIN STRUCTURE.') 808 FORMAT(' ',I3,' ANIMALS NO DATA: INFECTION AROUND ELECTRODE.') 809 FORMAT(' ',I3,' ANIMALS NO DATA: MISC!!!') 810 FORMAT(' ',I3,' ANIMALS PRODUCE DATA. TADA!') 811 FORMAT(' WITH BIPHASIC CURRENT CURRENT DIRECTION IGNORED.') 812 FORMAT(' WITH SINE WAVE CURRENT IS CONTINUOUS:IGNORE PULSE.') 813 FORMAT(' TEMPERATURE TOO LOW FOR SUBJECTS TO SURVIVE AT ALL.') 814 FORMAT(' ',I3,' ANIMALS DIED FROM A DRUG OVERDOSE.') 815 FORMAT(' ',I3,' ANIMALS DIED FROM THE LESION PROCEDURE.') 816 FORMAT(' ANESTHETIC DOSE TOO LOW TO OPERATE.') 817 FORMAT(/' MISSING DATA DUE TO EQUIPMENT OVERLOAD') 999 RETURN C C C C C 9999 ENTRY MINIT C OPEN(UNIT=IDEV3,ACCESS='SEQIN',FILE=MODNM2,DIRECTORY='4030,11') C READ(IDEV3,850)PRESS,SIGN,OVERD,DREFF,DARK,FIRSTT,FIRSTE 1,LMFB,EXPER2 850 FORMAT(24(10F5.1/),5F5.4,5F6.0/11(12F5.3/),8F5.3/ 12(12F5.2/),11I3/6F10.7/5F10.7/7(7F5.3/),I5) RETURN END