NY/*JOBPARM REGION=1024,LINES=19}// EXEC FORTGCLG //FORT.SYSIN DD *ChCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C0xPROGRAMIMSEGP.DET026.C0pVERSION3.414-SEP-820L6.1C0000T6.38C0000S7.\ C0000V7.001HCiDEVELOPEDFROMPROGRAMSISODATA.DETANDIMSEGP.COMMON.D7.01 C0000U7.02q C0000U7.03CcTHE"IMSEGP"PROGRAMSAREFORSEGMENTINGMULTI-CHANNELDIGITAL7.1CIMAGES.(FORSINGLE-CHANNELDATATHE"IMSEG1"PROGRAMSMAYBEc7.2&CUSED.)00007.3 CIMSEGP.COMMONUSESDISTANCEINTHEMETRICOFTHEESTIMATEDb9.+CCOMMONCOVARIANCEMATRIX.IMSEG.DETUSESDIFFERENTCOVARIANCE10.CMATRICES,WITHADJUSTMENTBYTHEDETERMINANTS, I.E.,ITUSESTHE11.0CESTIMATEDLOGLIKELIHOODFORTHEGAUSSIAN MODELWITHDIFFERENT12.QCCOVARIANCEMATRICES.00013.){CfMANUAL MODE:NUMBEROFCLASSESANDINITIAL MEANSARE INPUT16./yCUSEPROGRAMSIMSEGP.*.AUTOTOTRYA RANGEOFNUMBERSOFT17.݌CCLASSES,WITHAUTOMATICSETTINGOFINITIALMEANS.018.C0000L19._-C0000L20.CPROGRAMMEDBY:00021.CDR.STANLEYL.SCLOVE0<312/996-2681022.`CDEPARTMENTOFQUANTITATIVEMETHODS%312/996-2676023.'CCOLLEGEOFBUSINESSADMINISTRATION0024.YCUNIVERSITYOFILLINOISATCHICAGO0025.d}CBOX 4348,CHICAGO,IL 606800026.hC0000L27.ChRESEARCHSUPPORTEDINPARTBY:0028.hC0000L29.mCONRCONTRACTN00014-80-C-0408,TASKNR042-4430l30.qCAROCONTRACTDAAG29-82-K-015500t31. C0000L32.CRESTRICTIONS(CANBEMODIFIED):0033.uCNR,NUMBEROFROWSOFARRAY,ATMOST256;0L34.x8CNC,NUMBEROFCOLUMNSOFARRAY,ATMOST256;0t35.\CIP,NUMBEROFCHANNELS,ATMOST20;0035.1HCK,NUMBEROFCLASSES,ATMOST29;0036.CITER,MAXIMUMNUMBEROFITERATIONS,20.0337.qC0000L38.CSUBROUTINE(S)CALLED:00T39.ĝCMATEQ, WHICH CALLS MATDT00340.C&C0000L41.GC0000L42.+CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 42.1 C0000M42.2K0 C0000M42.3OQDIMENSIONX(56,56,10),SUM(29,4)0043.{ DIMENSIOND(29),ICLUS(56,56)00c44.y C0000M44.1SDIMENSION IOTA(56),JOTA(29)00u44.2W C0000M44.3-DIMENSIONTITLE(18)00045. DIMENSIONNG(29),XMEAN(29,4)00l46.2DIMENSIONFMT(18)00047.6DIMENSIONSS(29,4,4),SSD(29,4,4)0048.DIMENSION WGSS(4,4)00049.YDIMENSIONVARHAT(4,4),WGMS(4,4)0050.}DIMENSIONICLSOL(56,56)00L51. DIMENSION XMIN(4),XMAX(4)0052.CCHANGENO.OFCHANNELSFROM4TO200052.1 hCANDNO.OFROWSANDCOLUMNSFROM56TO2560M52.2׀CINDIMENSIONSTATEMENTSWHENALLOWABLEREGIONISINCREASED.D52.3 DIMENSIONIV(20,20)00053. DIMENSIONP(20,20)00054.C0000L55.DIMENSIONA(20,20)00056.8C0000L57.\DIMENSIONET(29)00058.HDIMENSIONPG(29,20,20)00059. C0000M59.1qDIMENSIONNT(29,29,29),IRSUM(29,29),TP(29,29,29)059.2DIMENSIONPROB(29)00059.3 C0000M59.4 &C0000L60.DOUBLEPRECISIONSS,SUM00361.DOUBLEPRECISIONWGSS,SSD00L62.<DOUBLEPRECISIONVARHAT00363.|0DOUBLEPRECISIONP00064.QDOUBLEPRECISIONDET00065.{DOUBLEPRECISIOND00066.$yDOUBLEPRECISION XMEAN00067.^ DOUBLEPRECISIONTEMPIV,TEMPJV00c68.DOUBLEPRECISIONF00069.-DOUBLEPRECISIONCF00070.,C0000L71.?DOUBLEPRECISIONA00072.C0000L73.DOUBLEPRECISIONET00074.:YDOUBLEPRECISIONPG00075.=}DOUBLEPRECISIONTP,PROB00U75.1bC0000L76.fCIVISAWORK ARRAYFORSUBROUTINEMATEQ.0L77.hC0000L78. C0000M78.1k+CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 78.2o C0000M78.3 C0000M78.4CCONTROLCARDS:000/C79.~8 C0000~C79.1v\CDATASET TITLE000C79.2zHCNUMBEROF ROWS,NR,INFORMAT(3X,I3)0C79.3[CNUMBEROFCOLUMNS,NC,INFORMAT(3X,I3)0C79.4qCNUMBEROFCHANNELS,IP,INFORMAT(3X,I2)0jC79.5 LCFMT,INFORMAT(18A4), E.G.,(1X,F4.1)0AC79.6Ce"FMT"WILLALSOBEUSEDFOROUTPUT: ALLOWAT LEASTONE BLANKC79.7]&CbATTHEBEGINNINGFORCARRIAGECONTROL.0jC79.8A CDATA,ONENUMBERATA TIME,INFORMATSPECIFIEDBYFMTC79.9ECdDATAISINDEXEDBYCHANNEL,ROW,ANDCOLUMN.0C80.ICfCOLUMNCHANGESFIRST,THENROW,THENCHANNEL.0C80.10 C0000~C80.2}QCK,NUMBEROFCLASSES,INFORMAT(2X,I2)0C80.3M{CKINITIALMEANS,INFORMATSPECIFIEDBYFMT0C80.4QyCgINITIAL MEANSAREINDEXEDBY CLASSANDCHANNEL.1C80.5CgCHANNELCHANGESFIRST,THENCLASS.0JC80.6\ C0000~C80.7U- C0000~C80.8Y+CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 80.9C0000L81.0 C0000M81.14READ(5,1050) TITLE00081.27Y C0000M81.3}CWRITEPROGRAMINFORMATION.00D81.4WRITE(6,24000)000d81.5WRITE(6,72500)000e81.51hWRITE(6,42000) TITLE00081.6 C0000M81.7CREADNUMBERSOFROWSANDCOLUMNS,NRANDNC.0m81.8 READ(5,12000)NR00081.9READ(5,12000)NC00082.WRITE(6,38000)NR00082.18WRITE(6,40000)NC00082.2\CREADNUMBEROFCHANNELS(I.E.,VARIABLES),IP.082.3HREAD(5,1105)IP00082.4WRITE(6,1110)IP00082.5q C0000M82.6CREADDATAFORMAT.00082.7READ(5,1050)FMT00082.8& C0000M82.9CREAD DATA.000L83.DO300 ICHAN=1,IP00U83.10DO300I=1,NR00M83.200DO300J=1,NC00d83.3(Q0READ(5,FMT)X(I,J,ICHAN)0D83.4{0IF(I.EQ.1 .AND.J.EQ.1)GOTO100m83.5y0GOTO20000U83.6 1000CONTINUE00483.7)0DO200IVAR=1,IP0483.8/-0XMAX(IVAR)=X(1,1,IVAR)0m83.9ݫ0XMIN(IVAR)=X(1,1,IVAR)0l84. 2000CONTINUE00u84.1_0DO300IVAR=1,IP0484.20IF(X(I,J,IVAR).LT.XMIN(IVAR))XMIN(IVAR)=84.3YX0X(I,J,IVAR)0084.4}0IF(X(I,J,IVAR).GT.XMAX(IVAR))XMAX(IVAR)=84.5@X0X(I,J,IVAR)0084.6 3000CONTINUE00u84.7dhWRITE(6,72010)000d84.8hWRITE(6,FMT)(XMIN(IVAR),IVAR=1,IP)0084.9WRITE(6,74010)000c85. WRITE(6,FMT)(XMAX(IVAR),IVAR=1,IP)0085.1mCREADK,NUMBEROFCLASSES.00u85.2qREAD(5,10000)K00085.38WRITE(6,26000)K00085.4\ C0000M85.5tHWRITE(6,20000)000d85.6xCREADINITIAL MEANS00085.7qDO400IG=1,K00085.80DO400IVAR=1,IP00u85.90READ(5,FMT)XMEAN(IG,IVAR)0386.& 400CONTINUE00086.1 C0000M86.2DO500IG=1,K00086.3CWRITE(6,22000)00086.4F0WRITE(6,FMT)(XMEAN(IG,IVAR),IVAR=1,IP)086.5Q 500CONTINUE000u86.6{ C0000M86.7KyCSETCONSTANTS.000m86.8O C0000M86.9/N= NR*NC000T87.-XN=N0000 87.01SDO600INTEG=1,NC00087.1WIOTA(INTEG)= INTEG00U87.2 600CONTINUE000u87.3DO705J=1,K00087.41YJOTA(J)=J00087.55} 705CONTINUE000u87.69 C0000M87.7PI=3.1415927000m87.8hC!PARAMETERSFOR MODELWITHCOMMONCOVARIANCEMATRIXAREM142. CKMEANVECTORSANDONECOVARIANCEMATRIX.0143.NOPARM=K*IP+IP*(IP+1)/200u144. C!PARAMETERSFOR MODELWITHDIFFERENTCOVARIANCEMATRICESAREd145.CKMEANVECTORSANDKCOVARIANCEMATRICES.0146.  NPRMDF=K*IP+K*IP*(IP+1)/200d147.8 C0000D148.\C^NON-DISTRIBUTIONALPARAMETERS:00148.01HCK**2-BY-KTRANSITION PROB.MATRIX00148.02CGIVES(K**2)(K-1)FREETRANSITIONPROBABILITIES0148.021qC0000G148.022NOPARM=NOPARM+(K**2)*(K-1)00148.023NPRMDF=NPRMDF+(K**2)*(K-1)00148.04& C0000F148.05 C0000E148.1CFOR FIRSTITERATION,MARGINALDISTRIBUTIONOFLABELS, PROB,e148.2 CIS TAKENTOBEUNIFORM.00V148.30 C0000E148.4QDO145IG=1,K000n148.5.{145PROB(IG)= 1.0/K000148.6|y C0000E148.7 C0000E148.8ITER=1000U149.$-C!ITERATIONS BEGINHERE000 149.1^ C0000E149.2 700CONTINUE000U150.IF (ITER.EQ.1)GOTO90000m151.,DO800I=1,NR000152.>YDO800J=1,NC000 152.1}ICLSOL(I,J)=ICLUS(I,J)00d153. 800CONTINUE000m154. 5:CSAVEPREVIOUS VALUEOF-2LOGMAXLIKELIHOOD:0155.=hXMN20L=XMN2LL000155.1b 900CONTINUE000U156.f C0000E156.1 CCOMMENCEDISTANCECOMPUTATIONS.00156.2 C0000E156.3kDO1700I=1,NR000157.n8DO1700J=1,NC000 157.1r\C!INITIALIZEDISTANCES(TOBEACCUMULATED)AT ZERO.0157.3H0DO1000L=1,K00D158.~0D(L)=0.0000159.vq1000CONTINUE000160.zCFOR FIRSTITERATION,EUCLIDEANDISTANCEISUSEDBECAUSED161.[CNOCOVARIANCEMATRIXISYETAVAILABLE. AFTERTHE FIRSTD162.&C!ITERATION,DISTANCEWILLBE TAKENINTHEMETRICOFTHEM163.CGROUPCOVARIANCEMATRIX.00M164. C0000D165.]IF (ITER.GT.1)GOTO120000166.{00DO1100 L=1,K00U167.DQ0DO1100IVAR=1,IP00168.H{0D(L)=D(L)+(XMEAN(L,IVAR)-X(I,J,IVAR))**2169.y11000CONTINUE00U170.}GOTO1500000171.M1200CONTINUE000m172.Q-0DO1400 L=1,K00U173.0DO1300IVAR=1,IP00174.\0TEMPIV=XMEAN(L,IVAR)-X(I,J,IVAR)0175.U0DO1300JV=1,IP00176.Y0TEMPJV=XMEAN(L,JV)-X(I,J,JV)4177.Y C0000D178.}0D(L)=D(L)+TEMPIV*PG(L,IVA R,JV)*TEMPJV179.3 C0000D180.713000CONTINUE00m181.h C0000E181.1 C0000E182.11400CONTINUE000183. 1500CONTINUE000m184. DO112L=1,K000e184.1CFOR FIRSTITERATION,CLASSIFICATIONISSIMPLYBY0184.28CMINIMUMDISTANCE.OTHERWISE,THETRANSITION0e184.3\C)PROBABILITIESANDDETERMINANTSENTER:05184.4HIF (ITER.EQ.1)GOTO11200n184.5CUPTONOW,D(L)IS(SQUARED)DISTANCE05184.6qCITMUSTNOWBEMODIFIEDTOPROBABILITY,0O184.61CFORMULTIPLICATIONBYTHETRANSITIONPROBABILITIES.0184.62CADDDETERMINANTOFCOVARIANCEMATRIXTODISTANCE:0184.621&0D(L)=D(L)+ ET(L)00g184.622ARG=-D(L)/2.0000184.63BIF(ARG.LT.-180.2)GOTO61505184.7BGOTO610000184.80615D(L)=0.0000E184.9QGOTO616000O184.91{610CONTINUE000W184.92y C0000F184.93(CMOVEFROMLOGPROB SCALETOPROBSCALE:0P184.931D(L)=EXP(ARG)000184.94-CD(L)ISNOLONGERADISTANCE:00g184.941CITIS(PROPORTIONALTO)THEBELONGINGPROBABILITY0184.942)616CONTINUE000W184.95/ C0000F184.96IF(I.EQ.1.OR.J.EQ.1)GOTO6050f184.97YIG1=ICLUS(I,J-1)000184.98%}IG2=ICLUS(I-1,J)000184.99׼D(L)=-TP(IG1,IG2,L)*D(L)00184.991 GOTO112000P184.992h605CONTINUE000X184.993@CCLASSIFYBOUNDARYOBSERVATIONS00184.994D(L)=-D(L)*PROB(L)000184.995d 112CONTINUE000X184.996hC0000G184.997F=D(1)000185.8ICLUS(I,J)=1000186.l\0DO1700L=2,K00D187.pH0IF(D(L)-F)1600,1700,17000188.1600F=D(L)000189.q0ICLUS(I,J)=L00M190.t1700CONTINUE000191.xWRITE(6,1003)000d192.& C0000D193.WRITE(6,1004)(IOTA(J),J=1,NC)00194.DO315I=1,NR000195.WRITE(6,1005)I,(ICLUS(I,J),J=1,NC)0N195.10 315CONTINUE000U196.QDO1800IG=1,K000 196.1B{NG(IG)=0000197.Fy0DO1800IVAR=1,IP00d198.ˌ0SUM(IG,IVAR)=0.000m199.0DO1800JV=1,IP00200.K-0>SS(IG,IVAR,JV)=0.000201.O0?SSD(IG,IVAR,JV)=0.000202.18000CONTINUE00U203.DO1900I=1,NR000204.SDO1900J=1,NC000 204.1VYIGROUP=ICLUS(I,J)00M205.Z}NG(IGROUP)=NG(IGROUP)+100206.0DO1900IVAR=1,IP00d207.10SUM(IGROUP,IVAR)=SUM(IGROUP,IVAR)+X(I,J,IVAR)d208.5h0DO1900JV=1,IP00209.9׫SS(IGROUP,IVAR,JV)=SS(IGROUP,IVAR,JV)+X(I,J,IVAR)*X(I,J,JV)210.19000CONTINUE00U212.  C0000E212.1 WRITE(6,1045)(NG(IG),IG=1,K)00212.3  C0000E212.4 8DO2000IVAR=1,IP00M213. \0DO2000JV=1,IP00u214. H0WGSS(IVAR,JV)=0.000d215. 2000CONTINUE000216. q C0000D217.  C0000D218. DO2300IG=1,K000219. &IF(NG(IG).EQ.0)GOTO21000U220. GOTO2200000221. 2100WRITE(6,78000)IG00U222. GOTO3500000223. 02200CONTINUE000m224. Q0DO2300IVAR=1,IP00d225. {CCOMPUTEMEANVECTORS000 225.1 Z y0XMEAN(IG,IVAR)=SUM(IG,IVAR)/NG(IG)0226. CCOMPUTESUM-OF-PRODUCTSMATRICES00226.1 .0DO2300JV=1,IP00227. |-0CF=SUM(IG,IVAR)*SUM(IG,JV)/NG(IG)U228. 0?SSD(IG,IVAR,JV)=SS(IG,IVAR,JV)-CFD229. 23000CONTINUE00U230. $ C0000D231. ^ C0000D232. Y C0000D233. }CPOOL:0000234. C0000D235. >DO2400IG=1,K000236. h0DO2400IVAR=1,IP00d237. 0DO2400JV=1,IP00238. :0_WGSS(IVAR,JV)=WGSS(IVAR,JV)+SSD(IG,IVAR,JV)239. = 24000CONTINUE00U240. b C0000D241. f C0000D242. i8CCOMPUTEVARHAT,MLEOFCOMMONCOVARIANCEMATRIX:0243. \DO2500IVAR=1,IP00M244. jH0DO2500JV=1,IP00u245. n0VARHAT(IVAR,JV)=WGSS(IVAR,JV)/N0m246. rq2500CONTINUE000247. CCOMPUTEDET(VARHAT)000251.1 ~IDET=1000U252. v&NRS1=0000D253. zDO2805IVAR=1,IP000253.1 [DO2805JV=1,IP000253.2  A(IVAR,JV)=VARHAT(IVAR,JV)00n253.3 02805CONTINUE000V253.4 QCALLMATEQ(A,IP,20,JFLG,DET,IDET,IV ,NRS1,P,20)0254. { C0000D255. {yCGENERALFORMOFCALLIS:00M256. DCCALLMATEQ(A,M,N,JFLG,DET,IDET,IV,NRS1,P,LL)0d257. HCSEESUBROUTINELISTINGFORFULLEREXPLANATION.0258. -C$DET(VARHAT)=DET*10.0**IDET00m259. }C(ACTUALDET.=DET*10**IDET)00n259.1 M C0000D260. QIF (JFLG.GT.0) WRITE(6,62000)JFLG04261. C0000D263. Y XIDET=IDET000u268. T}XLGDET=DLOG(DET)+XIDET*ALOG(10.0)04269. XXMN2LL=N*(IP*ALOG(2.0*PI)+IP+XLGDET)0u270. C0000D271. h C0000F272.01 3DO250IG=1,K000272.02 7_WRITE(6,1020)IG,(XMEAN(IG,IVAR) ,IVAR=1,IP)0272.03 250CONTINUE000W272.04 IF (ITER.EQ.1)GOTO60100f272.05DO555I=1,NR000272.068DO555J=1,NC000272.07\_IF(ICLUS(I,J).EQ.ICLSOL(I,J))GOTO5550o272.08HGOTO601000O272.09555CONTINUE000X272.091q GOTO530000P272.092601CONTINUE000X272.093C0000G272.094&CCOMPUTETRANSITIONPROBABILITYMATRIX:0V272.1 C0000E272.2DO3900I1=1,K000 272.30DO3900I2=1,K00272.400DO3900J=1,K00272.5Q0,NT(I1,I2,J)=000272.6{39000CONTINUE00V272.7yDO4000I=2,NR000 272.80DO4000J=2,NC00272.90IM1=I-1000273.-0JM1=J-1000273.1(0INORTH=ICLUS(IM1,J)00273.20 IWEST=ICLUS(I,JM1)00273.30IY=ICLUS(I,J)00E273.40NT(IWEST,INORTH,IY)=NT(IWEST,INORTH,IY)+1E273.5*Y4000CONTINUE000273.6-}DO4100I1=1,K000 273.70DO4100I2=1,K00273.80IRSUM(I1,I2)=000273.9%h4100CONTINUE000274.DO4200I1=1,K000 274.10DO4200I2=1,K00274.2 0DO4200J=1,K00274.3@0%IRSUM(I1,I2)=IRSUM(I1,I2)+NT(I1,I2,J)274.442000CONTINUE00V274.5c8DO4500I1=1,K000274.6g\0DO4500I2=1,K00N274.7H0XDENOM=IRSUM(I1,I2)00e274.8ŵ0IF(XDENOM.EQ.0.0)GOTO43000274.9lq0GOTO4400000275.p4300XDENOM=K000 275.14400CONTINUE000275.2& C0000E275.3t0DO4500J=1,K00275.4x0XNUM=NT(I1,I2,J)00275.5CIF THEREARENOTRANSITIONSFROM(I1,I2),THENTP(I1,I2,J)n275.60CISSET EQUALTO ZERO,FORALLJ=1,2,...,K.0e275.7Q C0000E275.8{45000TP(I1,I2,J)=XNUM/XDENOM0E275.9y4600CONTINUE000U276.، C0000E276.1BCCOMPUTEMARGINALDISTRIBUTIONOFLABELS:0E276.2F-DO4700IG=1,K000 276.3˫PROB(IG)=NG(IG)005276.4PROB(IG)=PROB(IG)/N00276.5K4700CONTINUE000n276.6O C0000E276.7RY TRANS=0.0000276.8}DO4800I1=1,K000276.90DO4800I2=1,K00M277.V0DO4800 J=1,K00e277.1Zh0ITEST=NT(I1,I2,J)00277.2 0IF(ITEST.EQ.0)GOTO48000277.310IF(TP(I1,I2,J).EQ.0.)GOTO48005277.45 0TRANS= TRANS+NT(I1,I2,J)*DLOG(TP(I1,I2,J))277.5948000CONTINUE00V277.6 C0000E277.78TRANS=-2.0*TRANS000277.8\C0000G277.801HCWRITETRANSITIONPROBABILITIES:00277.802 WRITE(6,1065)000g277.803qWRITE(6,1067)(JOTA(JAY),JAY=1,K)00277.804DO5000I1=1,K000277.8050DO5000I2=1,K00277.806 ŧ&0 WRITE(6,1066)I1,I2,(NT(I1,I2,J),J=1,K)7277.8075000CONTINUE000277.808WRITE(6,1070)000g277.809WRITE(6,1067)(JOTA(JAY),JAY=1,K)00277.810DO5100I1=1,K000277.811Q0DO5100I2=1,K00P277.812{0 WRITE(6,1078)I1,I2,(TP(I1,I2,J),J=1,K)7277.813y5100CONTINUE000277.814WRITE(6,88000)(PROB(IG),IG=1,K)00277.815WRITE(6,1201) TRANS000277.816-WRITE(6,30000)XMN2LL00X277.817 C0000F277.82.CACCOUNTFORLABELSOFBORDEROBSERVATIONS:0w277.83|FIRST=0.0000277.84DO516J=1,NC000278.01YLABEL1=ICLUS(1,J)000278.02!}PROBAB=PROB(LABEL1)000 278.1;FIRST= FIRST+ALOG(PROBAB)00o278.11516CONTINUE000W278.12hDO517I=2,NR000278.13LABEL1=ICLUS(I,1)000278.14>PROBAB=PROB(LABEL1)000278.15 FIRST= FIRST+ALOG(PROBAB)00o278.16517CONTINUE000W278.17:FIRST=-2.0*FIRST000278.2'8 C0000E278.3a\ C0000F278.31eHCCOMPUTE MODELSELECTIONCRITERIA:00278.32iCCOMPUTEVALUESCORRESPONDINGTONEWSEGMENTATIONANDOLD278.4qC!COVARIANCEMATRIX (OMITONFOR FIRSTITERATION)0278.401jIF (ITER.EQ.1)GOTO490000v278.6nCFOR MODELWITHCOMMONCOVARIANCEMATRIX0W278.61r&C(NOTOPTIMIZEDINTHISPROGRAM;HOWEVER,ITIS CLEAR6278.62CTHATONESHOULDUSETHECOMMON-COVARIANCE-MATRIX MODELIFx278.621~CTHEVALUESOFTHE MODELSELECTIONCRITERIAHEREFORTHAT MODEL278.63vCARELESSTHAN THOSEFORTHE MODELWITHDIFFERENTCOVARIANCEg278.631y0C!MATRICES):000G278.632Q C0000F278.65{AICOLD=XMN2OL+ TRANS+ FIRST+2.0*NOPARM0278.7ySCHOLD=XMN2OL+ TRANS+ FIRST+ALOG(XN)*NOPARM0278.71WRITE(6,1267)000f278.72WRITE(6,1260)AICOLD000 278.8{-WRITE(6,1261)SCHOLD000278.81D4900CONTINUE000V278.9HAIC=XMN2LL+ TRANS+ FIRST+2.0*NOPARM0279.WRITE(6,1265)AIC000279.1}SCH=XMN2LL+ TRANS+ FIRST+ALOG(XN)*NOPARM0279.11LYWRITE(6,1266)SCH000279.12P} C0000E279.2 C0000F280.81CCOMPUTEPRECISIONMATRICESFORALLGROUPS.0w280.82Th C0000F280.83XTERM=0.0000W280.84DO3000 L=1,K000280.841 0DO2800IVAR=1,IP00x280.84230DO2800JV=1,IP00280.84370A(IVAR,JV)=SSD(L,IVAR,JV)/ NG(L)0280.844828000CONTINUE00X280.845\IDET=1000p280.846HNRS1=0000p280.847CALLMATEQ(A,IP,20,JFLG,DET,IDET,IV ,NRS1,P,20)0280.848q0DO2900IVAR=1,IP00x280.8490DO2900JV=1,IP00280.85 0_PG(L,IVAR,JV)=P(IVAR,JV)0280.851&29000CONTINUE00X280.852C0000G280.853ET(L)=DLOG(DET*10**IDET)00280.854TERM=TERM+NG(L)*ET(L)00280.855 03000CONTINUE000p280.856QCCOMPUTE MODELSELECTIONCRITERIA BASEDONDIFFERENTCOVARIANCE280.857{CMATRICES:000P280.858yC0000G280.859C$PARAMETERS:000280.9CKMEANVECTORSOFDIMENSIONPANDKP-BY-PCOVARIANCEMATRICES,280.91-CWHEREPISTHENUMBEROFCHANNELS00280.92XM2LLD=N*IP*ALOG(2*PI)+N*IP+TERM0W280.93XM2LLD=XM2LLD+ TRANS+ FIRST00280.94SCHD=XM2LLD+ALOG(XN)*NPRMDF00280.95(AICD=XM2LLD+2.0*NPRMDF00280.96&YWRITE(6,31900)000280.961}WRITE(6,32000)XM2LLD00W280.97WRITE(6,76000)NPRMDF006280.98*WRITE(6,74000)AICD000280.99-hWRITE(6,74005)SCHD000281. C0000F281.01 C0000F281.02% WRITE(6,34000)ITER000281.03DO3100IG=1,K000281.04WRITE(6,28000)IG,(XMEAN(IG,IVAR) ,IVAR=1,IP)6281.0583100CONTINUE000o281.06#\ITER=ITER+1000281.07"HIF (ITER.GE.21)GOTO520000281.08cCUNLESS20ITERATIONS,HAVEALREADYBEENPERFORMED,GOBACKo281.09gqCANDDOANOTHER.000281.1GOTO700000O281.11ŝ5200 WRITE(6,1160)000f281.12l& 5300STOP0000281.13p C0000F281.14 C0000F281.15530CONTINUE000W281.16s0 C0000F281.17wQCOUTPUTTOBEWRITTENUPONCONVERGENCE:0W281.18{ C0000F281.19yWRITE(6,1120)ITER000281.2WRITE(6,1077)000f281.22DO5500IVAR=1,IP00N282.1-0DO5500JV=1,IP00v282.2 ث0WGMS(IVAR,JV)=WGSS(IVAR,JV)/(N-K)0282.3B5500CONTINUE000282.4FDO5600IVAR=1,IP005282.5WRITE(6,1080)(WGMS(IVAR,JV),JV=1,IP)0e282.6Y5600CONTINUE000n282.7J} C0000F282.71NWRITE(6,86005)000282.711RDO5510L=1,K000282.72hDO5505IVAR=1,IP000282.73DO5505JV =1,IP000282.74V_WGMS(IVAR,JV)=SSD(L,IVAR,JV)/(NG(L)-1)0F282.75Z 5505CONTINUE000W282.76>WRITE(6,86000)L000282.771DO5510IVAR=1,IP000282.7848_WRITE(6,1080)(WGMS(IVAR,JV),JV=1,IP)06282.798\5510CONTINUE000X282.791HC0000G282.792 C0000E282.8q C0000E282.9 C0000E283.8 C0000E283.9& C0000D358.WRITE(6,1220)000e358.1 C0000D359.3500STOP0000360.0 C0000E360.1Q1003FORMAT(//1X,'SEGMENTATION:'/)00e360.2{1004FORMAT(/,1X,'ROW:gCOLUMN:'/,4X,(40I3/))V360.3y1005FORMAT(1X,I3,(40I3/))00N360.41065FORMAT(/1X,'TRA NSITIONS'/)00360.4011110FORMAT(/1X,'NUMBEROFCHANNELS=',I2/)0O360.41- C0000E360.5 10000FORMAT(2X,I2)000m361. 12000FORMAT(3X,I3)000m362. 20000FORMAT(//1X,'IN ITIALMEANS'/)00d366. 22000FORMAT(1X,(8E13.5/)//)00U367.Y(24000 FORMAT('1','###################################################', 368.+}^X//,1X,'PROGRAMIMSEG.DET'/00369.X,1X,'FOR IMAGESEGMENTATION'/00370.ȟ$X,1X,'USINGDISTANCEINTHEMETRICSOFTHECOVARIANCEMATRICES'/371.!h;X,1X,'ADJUSTEDBYTHEDETERMINANTS0K'/372.;צ^X//,1X,'DEVELOPEDANDPROGRAMMEDBYDR.STANLEYL.SCLOVE'u373.ߨ^X//,1X,'VERSION3.414-SEP-82'//)00374.  26000FORMAT('1',//,1X,'K=',I2,'CLASSES')0U375. 28000FORMAT(1X,'MEANVECTORFOR CLASS',I2,':',(8E13.5/))M376.> 30000FORMAT(/,1X,'MINUS2LOGLIKELIHOODFOR MODELWITHCOMMON',d377.8X'COVARIANCEMATRIX=',E13.5//)00378.\+31900 FORMAT(//1X,'FOR MODEL BASED ON DIFFERENT COVARIANCE MATRICES:'//) 378.1`H< 32000FORMAT(/,1X,'MINUS2LOGLIKELIHOODFOR MODELWITHDIFFERENT',379.'*X'COVARIANCEMATRICES=',E13.5//)00380.aq 34000FORMAT(///,1X,'ITERATION',I2,//)00381.e1020FORMAT(1X,'MEANVECTORFOR CLASS',I2,':',(8E13.5/))O381.01i1045FORMAT(/,1X,'NUMBERS:',(9I13/)/)00381.05&1050FORMAT(18A4)000v381.1j1105FORMAT(3X,I2)000n381.2n<1077FORMAT(///,1X,'COMMONCOVARIANCEMATRIX(DIVISORISDF):',//)381.201r 1080FORMAT(1X,(8E13 .5/))000381.210.1201FORMAT(//,'CONTRIBUTIONOFTRANS.PROBS.TOLOGLIKELIHOOD=',381.3QXE15.5/)0005381.4u{ 38000FORMAT(1X,'NUMBEROFROWS= ',I4)00383.yy 40000FORMAT(1X,'NUMBEROFCOLUMNS=',I3/)04384.ڌ 42000FORMAT(1X,18A4)000385. 46005FORMAT(1X,2I4,(9I7/))00V387.1- 44005FORMAT(/1X,'TRA NSITIONS'/)00387.2 48005FORMAT(/9X,(9I7/))000387.3(48000 FORMAT(///,1X,'COMMON COVARIANCE MATRIX (UNBIASED ESTIMATE):',//) 388.{ 50000FORMAT(1X,8E13.5/)000389.D 52000FORMAT(1X,I4,1X,I2)000390.GY1067FORMAT(/9X,(9I7/))000390.001}1066FORMAT(1X,2I4,(9I7/))00W390.011070FORMAT(/1X,'TRANSITIONPROBABILITIES'/)0N390.1L1078FORMAT(1X,2I4,3 X,(9F7.4/))00391.Ph1120FORMAT(/1X,'CONVERGENCE:NOCASECHANGEDCLASSES AFTER',m394.׫$X'ITERATION',I2,'.SOMEADDITIONALRESULTSAREPRINTEDBELOW.'//)395. 62000FORMAT(/,1X,'JFLG=',I2,'.IFJFLG=0,COMPUTATIONOF DET',396.T X'WENT WELL;OTHERWISE, THEREWASTROUBLEORMATRIXWAS',m397.X^X'ILL-CONDITIONED.'//)000398.1160FORMAT(1X,'ROUTINEHASNOTCONVERGEDIN20ITERATIONS.STOP')403.81220FORMAT(/,1X,'PR OGRAM ENDEDSUCCESSFULLY.')0v403.12\(72000 FORMAT(1X,'AIC FOR MODEL WITH COMMON COVARIANCE MATRIX = ',F16.5/) 404.6H 72500FORMAT(/1X,'MODELWITHCOMMONCOVARIANCEMATRIXIS NOT'/404.02;X1X,'OPTIMIZEDINTHISPROGRAM;HOWEVER,ITISCLEAR'/X404.021qX1X,'THATONESHOULDUSEITIFHERETHEVALUESOFTHEMODEL-'/404.022;X1X,'SELECTIONCRITERIAFORTHAT MODELARELESSTHAN THOSE'/404.023X1X,'FORTHE MODELWITHDIFFERENTCOVARIANCEMATRICES,'/G404.024 T&!X1X,'WHICHISOPTIMIZEDHERE.'//)00404.025 72010FORMAT(/1X,'MINIMUMFOREACHCHANNEL:',/)0w404.07 74010FORMAT(/1X,'MAXIMUMFOREACHCHANNEL:',/)0w404.08| 72005FORMAT(1X,'SCHWARZCRITERIONFOR MODELWITHCOMMONCOVARIANCE',404.090X'MATRIX=',F16.5/)000404.1Q+ 74000FORMAT(1X,'AICFOR MODELWITHDIFFERENTCOVARIANCEMATRICES=',404.11{XF16.5/)0004406. y 74005FORMAT(1X,'SCHWARZCRITERIONFOR MODELWITHDIFFERENT',406.1*X'COVARIANCEMATRICES=',F16.5/)06406.111260FORMAT(/1X,'AIC BASEDONNEWLABELSANDOLDDISTRIBUTIONAL',406.45-)X'PARAMETERS:',E15.5)00W406.461261FORMAT(/1X,'SCHWARZCRITERION BASEDONNEWLABELSANDOLD',406.47^X'DISTRIBUTIONALPARAMETERS:',E15.5)0W406.481265FORMAT(/1X,'AIC BASEDONNEWLABELSANDNEWDISTRIBUTIONAL',406.5)X'PARAMETERS:',E15.5)005406.6Y1266FORMAT(/1X,'SCHWARZCRITERION BASEDONNEWLABELSANDNEW',406.7}^X'DISTRIBUTIONALPARAMETERS:',E15.5)05406.8<1267FORMAT(//1X,'FOR MODEL BASEDONCOMMONCOVARIANCEMATRIX',f406.81&X'(NOTOPTIMIZEDINTHISPROGRAM):'//)0O406.82h 76000FORMAT(/,1X,'NUMBEROFPARAMETERS=',I4//)0m407. 78000FORMAT(1X,'NOOBSERVATIONSIN GROUP',I3,'.STOP')0408.* 80000FORMAT(//,1X,'PROGRAM ENDEDNORMALLY.')0M409.-  86000FORMAT(//1X,'COVARIANCEMATRIXFOR CLASS',I2/)0412. 86005FORMAT(///1X,'CLASSCOVARIANCEMATRICES',0w412.01X'(DIVISORSAREONELESSTHANNUMBERINGROUP):'/)0412.02,8 88000FORMAT(/1X,'MAR GINALPROB.VECTOR:',(9F11.4/))0412.1?\ END0000413.HSUBROUTINEMATEQ(A,M,N,JFLG,DET,IDET,IV,NRS1,P,LL)0414.ֵC!SUBROUTINE MATEQISDMATEQFROMTHEUICCSUBROUTINELIBRARY.415.#q C0000D416."C0SUBROUTINEDMATEQ0D417.cC0*****************0D418.g&<CTHISROUTINEWILL SOLVEAREAL*8SYSTEMOFLINEAREQUATIONS,COMPUTE419..CTHEDETERMINANT,WITHOUTUNDERFLOWOROVERFLOW,OFAREAL*8MATRIX,420.CAND/ORINVERTAREAL*8MATRIX.00u421.lCCALLINGSEQUENCE:000422.o0CCALLDMATEQ(A,N,IA,JFLG,DET,IDET,IV,NRS,P,IP)WHERE;D423.Q CjA(INPUT)-ISTHEREAL*8MATRIXON WHICHTHEROUTINEIS424.{C0TO WORK.INTHEPROCESSOFCOMPUTATIONTHE425.syC0CONTENTSOFTHISMATRIXAREDESTROYED.d426.w CjN(INPUT)-ISANINTEGER*4VARIABLE WHICHSPECIFIESTHE427.C0ORDEROFTHEAMATRIX.0m428.-CkIA(INPUT)-ISANINTEGER*4VARIABLE WHICHSPECIFIESTHE429.C0ACTUALROWDIMENSIONOFAASDIMENSIONEDIN430.C0THECALLINGPROGRAM.IAMUSTBEGREATERTHAN431.C0OR EQUALTON.00432.CmJFLG(OUTPUT)-ISANINTEGER*4RETURNCODEVARIABLE.UPON433.AYC0RETURNFROMDMATEQIF;0m434.E}C0JFLG=0,ALLWENT WELL.0m435.IC0JFLG=1,THEAMATRIXWASSINGULARORNEAR436.C0SINGULARANDTHECOMPUTATIONS COULDNOTBE437.JhC0COMPLETED.THECONTENTSOFTHEVARIABLES438.NC0A,DET,IDETANDPAREMEANINGLESS.439.RClDET(OUTPUT)-ISAREAL*8VARIABLE WHICHCONTAINSTHE440. C0DETERMINANTOFA.(SEE IDET)0441.CmIDET(INPUT)-ISANINTEGER*4VARIABLE.ON INPUTIF;d442.VC0IDET=0,NODETERMINANTISCALCULATED.m443.Y8C0IDETNOT0,THEDETERMINANTOFAISCOMPUTED.444.\C0KONOUTPUTIDETCONTAINSTHE POWEROF10445.0HC0THATDETSHOULDBEMULTIPLIEDBYTOGIVETHE446.4C0CORRECT VALUEOFTHEDETERMINANT.I.E.d447.8qC0DET(A)=DET*10.0D0**IDET.0448.C0KIFDET(A)CANBECOMPUTEDWITHOUT UNDEROR449. C0OVERFLOW,THENIDET=0OTHERWISEIDETISSET450. &C0TOTHEPROPER VALUESOTHATNO UNDEROR OVER-451. C0FLOWWILL OCCURINCOMPUTINGDET.M452. CkIV(INPUT)-ISANINTEGER*4WORK ARRAY WHICHSHOULDBE453.  C0DIMENSIONEDAT LEASTIV(N).0454. 0 C0000D455. QClNRS(INPUT)-ISANINTEGER*4VARIABLEWITHTHEFOLLOWING456. {C0INTERPRETATION:00457. yC0NRS>0, SOLVEASYSTEMOFLINEAREQUATIONS458. C0WITHNRS RIGHTHANDSIDES.0459.  C0NRS=0,INVERTTHEAMATRIX.0460. S -C0NRS<0,ONLYCOMPUTETHEDETERMINANTOFA.461. C0INTHISCASEIDETMUSTBEDIFFERENTFROM0.462. CjP(INPUT)-ISAREAL*8 ARRAYWITHTHEFOLLOWINGINTER-463.  C0PRETATION:00464. C0IFNRS>0,THENPCONTAINSTHENRS RIGHTHAND465. YC0SIDESSTOREDBYCOLUMNS.INTHISCASEPMUST466. }C0BEDIMENSIONEDAT LEASTP(N,NRS).ONRETURN467. C0THECOLUMNSOFPAREREPLACEDBYTHERESPEC-468. C0TIVESOLUTIONS.00469. +hC0IFNRS=0,THENPMUSTBEDIMENSIONEDAT LEAST470. C0P(N,N).ONRETURNPWILLCONTAINTHEINVERSE471. C0OFA.00472. ! C0IFNRS<0,THENPNEEDONLYBEA DUMMYVARIABLE473. ;C0INTHISCASEPIS NEVERACCESSEDBYDMATEQ.474.  CkIP(INPUT)-ISANINTEGER*4VARIABLE WHICHCONTAINSTHE475. 8C0ACTUALROWDIMENSIONOFPASDIMENSIONEDIN476. \C0THECALLINGPROGRAM.IPMUSTBEGREATERTHAN477. _HC0OR EQUALTON.00478. CNOTE:IMMEDIATELYONRETURNFROMDMATEQTHECONDITIONCODE FLAG,479. qCvJFLG,SHOULDBEINTERROGATED.IFJFLG=1,THENTHEROUTINE480. `CvCOULDNOTCOMPUTEASOLUTION.00481. 'CMETHOD-THEALGORITHMUSEDISGAUSSIANELIMINATIONWITHPARTIAL482. a&C0000-1d483. eCrPIVOTING.INESSENCETHEROUTINEGENERATESAMATRIXLSUCH484. i C2-1000u485. CmTHATL*A=U, WHEREUISAN UPPERTRIANGULARMATRIX.THENIT486. 0CoSOLVESTHESYSTEMA*X=PBY MEANSOFTHEEQUIVALENTSYSTEM487. mQC0-1b-100U488. q{ClU*X=L*A*X=L*PBYBACKSUBSTITUTION.0489. y C2-1000u490. ȌCtTHELMATRIXCANBEWRITTENASAPRODUCTOFTHEFORMd491. uCk-10000492. y-CaL=L*P*....*L*P WHEREEACHPISAPERMUTATION493. ګC3N-1N-1a11JK0d494.  CfMATRIXOBTAINEDBYINTERCHANGINGATMOSTTWOROWSOFTHEd495. ChIDENTITYMATRIX.(THISREPRESENTSTHEINTERCHANGINGOFTWO496. CfROWS).THELMATRICESAREELIMINATIONMATRICES WHICHARE497. Y C0K000498. }CfCHOSENTOINTRODUCE ZEROSINTHELASTN-KENTRIESOFTHEK-TH499. CCfCOLUMNOFTHEMATRIX.00M500. G C0000D501. hC02-1b-104502. ClTHECALCULATIONSOFL*AANDL*PAREDONEBYPERFORMING503. LCTHEPERMUTATIONSONAANDPRESPECTIVELY.THEACTUALLANDP504. P C0000K/K505. CARENOTCOMPUTED.000506. CSUBROUTINESCALLED:DMATDT00U507. S8 CREFERENCE:0004508. W\CG.W.STEWART,INTRODUCTIONTOMATRIXCOMPUTATIONS,0509. HCACADEMICPRESS, 1973.000510. REAL*8A(N,1),DET,P(LL,1)00D511. 2qREAL*8DNORM,DEN,DMULT ,DSUM,DISIGN00512. 6DIMENSION IV(1)000513. NRS=NRS1000U514. &IF(NRS.EQ.0)IDET=1000515. >DISIGN=1.0D+00000d516. ,DET=0.0D+00000517. JFLG=00000518. 0 C0000D519. Q(C JFLG IS A TROUBLE FLAG.UPON EXIT IF JFLG=0 THEN THE MATRIX WAS PROCESS 520. {CWITHOUTTROUBLE.IFJFLG=1EITHERTHEMATRIXISSINGULARORTROUBLE521. y(C OCCURED.ISIGN=-ISIGN EVERY TIME A ROW IS INTERCHANGED.THIS IS USED TO 522. CINSURETHATTHEDETERMINANTHASTHEPROPER SIGN.0d523.  C0000D524. -M1=M-10000525.  DO100 I=1,M000526. 100IV(I)=I000u527. IF (NRS)500,200,500000528. 200DO300 I=1,M000529. Y 0DO300 J=1,M004530. }300P(I,J)=0.0D+0000M531. DO400 I=1,M000532. 400>P(I,I)=1.0D+00000533. hNRS=M0000534. < C0000D535. &(C INSTEAD OF ACTUALLY INTERCHANGING ROWS A POINTER ARRAY IS USED TO KEEP 536. C TRACKOFTHEROWPOSITIONS.00U537.  C0000D538. *C BEGINELIMINATION LOOP.000539. ^8 C0000D540. \500DO1200K=1,M1000541. HICOL=K000542. ,IPCOL=K000u543. ?q C0000D544. CSEARCHINGFORLARGESTELEMENTINABSOLUTE VALUEINCOLUMNK.545. ֝ C0000D546. #&DNORM=A(IV(K),K)000547. "IFLG=0000548. cKK=K+1000549. g0DO600J=KK,M00U550. 00IF(DABS(A(IV(J),K)).LE.DABS(DNORM))GOTO600551. Q0IFLG=1000552. k{ 0IPCOL=IV(J)000553. oy0DNORM=A(IPCOL,K)00554. 600CONTINUE000555. C0000D556. s-.CIFIFLG=0NOROWINTERCHANGETOOKPLACE.IFIFLG=1AROWINTERCHANGE557. ŧ wCTOOK PLACEANDTHEPOINTER ARRAYIVMUSTBEUPDATED.0558.  C0000D559. IF(IFLG.EQ.0)GOTO80000d560. ISAVE=IV(ICOL)000561. YIV(ICOL)=IPCOL000562. } ICOL1=ICOL+1000563. ]0DO700L=ICOL1,M00564. A0IF(IV(L).EQ.IPCOL)IV(L)=ISAVE0565. Eh 700CONTINUE000566. IDISIGN=-DISIGN000567. 800IF(DNORM.EQ.0.0D+00)GOTO190004568. J C0000D569. NC BEGINELIMINATIONOFROW BELOWIV(K).DENISTHE PIVOTELEMENT.m570. R C0000D571. 8K1=K+1000572. \\0DO1100IM=K1,M00D573. UH C0000D574. Y(CBEFOREACTUALLYELIMINTINGWE CHECKTOSEEIFA(IV(IM),K)HASALREADY575. qCBEENANIHALATED.000u576. 0 C0000D577. 40IF(A(IV(IM),K).EQ.0.0D+00)GOTO11000578. 8& C0000D579. CCACULATEELIMINATIONFACTOR.00M580.  C0000D581. 0DMULT=-A(IV(IM),K)00m582. 0 C0000D583. Q CWENOWCALCULATE VALUEOF OTHERELEMENTSINROWIV(IM).0584. { C0000D585. y0DO900NN=K1,M00d586. 9000A(IV(IM),NN)=(DMULT*A(IV(K),NN))/DNORM+A(IV(IM),NN)587. 0IF(NRS.LE.0)GOTO110000588. -0DO1000IN=1,NRS00589. 10000P(IV(IM),IN)=(DMULT*P(IV(K),IN))/DNORM+P(IV(IM),IN)590. 1100CONTINUE000591. 1200CONTINUE000m592.  C0000D593. YCCALCULATE VALUEOFDETERMINANT.00u594. } C0000D595. IF(A(IV(M),M).EQ.0.0D0)GOTO190000596. DET=DISIGN000D597. hIF(IDET.NE.0)CALLDMATDT(A,N,M,DET,IV,IDET)0598. IF(DET.EQ.0.0D+00)GOTO190000599. IF(NRS.LE.0)GOTO200000M600. + C0000D601. (C WE START SOLVING RIGHT HAND SIDES.THE SOLUTION REPLACES THE RIGHT HAND 602. CVECTOR.0000603. 8 C0000D604. )\ 1300N1=M-10000605. /HDO1600JJ=1,NRS000606. ݵ C0000D607. qC BEGINBACKSUBSTITUTION.000608. _ C0000D609. P(IV(M),JJ)=P(IV(M),JJ)/A(IV(M),M)0M610. &0DO1500I=1,N100D611. ` 0DSUM=0.0D+00004612. '0DO1400 J=1,I00m613. a14000DSUM=DSUM-A(IV(M-I),M-J+1)*P(I V(M-J+1),JJ)u614. d01500P(IV(M-I),JJ)=(P(IV(M-I),JJ)+DSUM)/A(IV(M-I),M-I)d615. hQ1600CONTINUE000m616. {DO1800JJ=1,NRS000617. y0DO1700IND=1,M00D618. m1700A(IND,1)=P(IV(IND),JJ)00619. q0DO1800IND=1,M00D620. -1800P(IND,JJ)=A(IND,1)00m621. ȫRETURN0000622. u 1900JFLG=10000623. yIDET=00000624. 2000RETURN0000625. Y END0000626. }SUBROUTINEMATDT(A,IA,N,DET,IV,IDET)00627. C!SUBROUTINE MATDTISDMATDTFROMTHEUICCSUBROUTINELIBRARY.628. REAL*8A(IA,1),DET,B,LOG1600629. hINTEGER*4IV(1),K000630. C,EQUIVALENCE (B,K)000631. G %NUM=16777216000u632. ?LOG16=.120411998265592457D+0100d633. IF(A(IV(N),N).EQ.0.0D+00)GOTO30000634. L L=00000635. O8 DO100 I=1,N000636. \B=DABS(A(IV(I),I))00U637. HK=K/NUM-64000638. SL=L+K000D639. Wq100?DET=DET*(A(IV(I),I)/16.0D+00**K)04640. ,B=DABS(DET)000641. K=K/NUM-64000D642. 2&IW=L+K0000643. 6IF((IW.LT.-64).OR .(IW.GT.63))GOTO2000D644. ?DET=DET*16.0D+00**L000645. IDET=00000646.0GOTO400000M647. Q200DET=DET*16.0D+00**(-K)004648.{IDET=L+K000U649. y %B=IDET*LOG16000u650.IDET=B0000651. ?B=B-DFLOAT(IDET)000652.-?DET=DET*1.0D+01**B000653.GOTO400000M654.300DET=0.0D+00000655. IDET=00000656.400RETURN0000657.Y END0000658.}//GO.SYSIN DD *FISHER IRIS DATANR=015hNC=010 IP=04(1X,F3.1) 5.1VARIABLE1A1001< 4.9VARIABLE1A1002& 4.7VARIABLE1A10038 4.6VARIABLE1A1004\ 5.0VARIABLE1A1005$H 5.4VARIABLE1A1006^ 4.6VARIABLE1A1007q 5.0VARIABLE1A1008 4.4VARIABLE1A1009, 4.9VARIABLE1A1010?& 5.4VARIABLE1A1011 4.8VARIABLE1A1012 4.8VARIABLE1A1013# 4.3VARIABLE1A1014=0 5.8VARIABLE1A1015bQ 5.7VARIABLE1A1016f{ 5.4VARIABLE1A1017y 5.1VARIABLE1A1018 5.7VARIABLE1A1019k 5.1VARIABLE1A1020o- 5.4VARIABLE1A1021 5.1VARIABLE1A1022 4.6VARIABLE1A1023s 5.1VARIABLE1A1024w 4.8VARIABLE1A1025 4zY 5.0VARIABLE1A1026[} 5.0VARIABLE1A1027 5.2VARIABLE1A1028 5.2VARIABLE1A1029h 4.7VARIABLE1A1030] 4.8VARIABLE1A1031A 5.4VARIABLE1A1032E 5.2VARIABLE1A1033I 5.5VARIABLE1A1034 4.9VARIABLE1A1035}8 5.0VARIABLE1A1036M\ 5.5VARIABLE1A1037QH 4.9VARIABLE1A1038 4.4VARIABLE1A1039\q 5.1VARIABLE1A1040U 5.0VARIABLE1A1041Y 4.5VARIABLE1A1042& 4.4VARIABLE1A10430 5.0VARIABLE1A10444 5.1VARIABLE1A10458 4.8VARIABLE1A10460 5.1VARIABLE1A1047Q 4.6VARIABLE1A1048{ 5.3VARIABLE1A1049y 5.0VARIABLE1A1050 7.0VARIABLE1A2051 6.4VARIABLE1A2052- 6.9VARIABLE1A2053 5.5VARIABLE1A2054 6.5VARIABLE1A2055 5.7VARIABLE1A2056 6.3VARIABLE1A2057Y 4.9VARIABLE1A2058} 6.6VARIABLE1A2059 5.2VARIABLE1A2060 5.0VARIABLE1A2061h 5.9VARIABLE1A2062 6.0VARIABLE1A2063 6.1VARIABLE1A2064 5.6VARIABLE1A2065 6.7VARIABLE1A2066 5.6VARIABLE1A2067(8 5.8VARIABLE1A2068\ 6.2VARIABLE1A2069H 5.6VARIABLE1A2070 5.9VARIABLE1A2071)q 6.1VARIABLE1A2072/ 6.3VARIABLE1A2073ݝ 6.1VARIABLE1A2074& 6.4VARIABLE1A2075_ 6.6VARIABLE1A2076 6.8VARIABLE1A2077 6.7VARIABLE1A20780 6.0VARIABLE1A2079@Q 5.7VARIABLE1A2080{ 5.5VARIABLE1A2081dy 5.5VARIABLE1A2082h 5.8VARIABLE1A2083 6.0VARIABLE1A2084- 5.4VARIABLE1A2085m 6.0VARIABLE1A2086q 6.7VARIABLE1A2087 6.3VARIABLE1A2088 5.6VARIABLE1A2089tY 5.5VARIABLE1A2090x} 5.5VARIABLE1A2091 6.1VARIABLE1A2092 5.8VARIABLE1A2093h 5.0VARIABLE1A2094 5.6VARIABLE1A2095٨ 5.7VARIABLE1A2096 5.7VARIABLE1A2097C 6.2VARIABLE1A2098G 5.1VARIABLE1A20998 5.7VARIABLE1A2100\ 6.3VARIABLE1A3101KH 5.8VARIABLE1A3102O 7.1VARIABLE1A3103q 6.3VARIABLE1A3104 6.5VARIABLE1A3105S 7.6VARIABLE1A3106W& 4.9VARIABLE1A3107 7.3VARIABLE1A3108 6.7VARIABLE1A31092 7.2VARIABLE1A311050 6.5VARIABLE1A31119Q 6.4VARIABLE1A3112{ 6.8VARIABLE1A3113y 5.7VARIABLE1A3114 5.8VARIABLE1A3115 6.4VARIABLE1A3116 - 6.5VARIABLE1A3117 7.7VARIABLE1A3118  7.7VARIABLE1A3119 6.0VARIABLE1A3120 6.9VARIABLE1A3121Y 5.6VARIABLE1A3122} 7.7VARIABLE1A3123 6.3VARIABLE1A3124 6.7VARIABLE1A3125h 7.2VARIABLE1A3126 6.2VARIABLE1A3127 6.1VARIABLE1A3128 6.4VARIABLE1A3129  7.2VARIABLE1A3130 7.4VARIABLE1A31318 7.9VARIABLE1A3132.\ 6.4VARIABLE1A3133|H 6.3VARIABLE1A3134 6.1VARIABLE1A3135q 7.7VARIABLE1A3136$ 6.3VARIABLE1A3137^ 6.4VARIABLE1A3138& 6.0VARIABLE1A3139 6.9VARIABLE1A3140, 6.7VARIABLE1A3141? 6.9VARIABLE1A31420 5.8VARIABLE1A3143Q 6.8VARIABLE1A3144:{ 6.7VARIABLE1A3145=y 6.7VARIABLE1A3146b 6.3VARIABLE1A3147f 6.5VARIABLE1A3148- 6.2VARIABLE1A3149 5.9VARIABLE1A3150k 3.5VARIABLE21001o 3.0VARIABLE21002 3.2VARIABLE21003Y 3.1VARIABLE21004~} 3.6VARIABLE21005v 3.9VARIABLE21006z 3.4VARIABLE21007[h 3.4VARIABLE21008 2.9VARIABLE21009 3.1VARIABLE21010 3.7VARIABLE21011] 3.4VARIABLE21012A 3.0VARIABLE21013D8 3.0VARIABLE21014H\ 4.0VARIABLE21015H 4.4VARIABLE21016} 3.9VARIABLE21017Mq 3.5VARIABLE21018Q 3.8VARIABLE21019 3.8VARIABLE21020\& 3.4VARIABLE21021U 3.7VARIABLE21022Y 3.6VARIABLE21023 3.3VARIABLE210240 3.4VARIABLE210253Q 3.0VARIABLE210267{ 3.4VARIABLE21027y 3.5VARIABLE21028 3.4VARIABLE21029 3.2VARIABLE21030- 3.1VARIABLE21031 3.4VARIABLE21032 4 4.1VARIABLE21033 4.2VARIABLE21034 3.1VARIABLE21035Y 3.2VARIABLE21036} 3.5VARIABLE21037 3.6VARIABLE21038 3.0VARIABLE21039h 3.4VARIABLE21040 3.5VARIABLE21041 2.3VARIABLE21042 3.2VARIABLE21043 3.5VARIABLE21044 3.8VARIABLE210458 3.0VARIABLE21046\ 3.8VARIABLE21047H 3.2VARIABLE21048( 3.7VARIABLE21049q 3.3VARIABLE21050 3.2VARIABLE22051 3.2VARIABLE22052)& 3.1VARIABLE22053/ 2.3VARIABLE22054 2.8VARIABLE22055 2.8VARIABLE22056%0 3.3VARIABLE22057Q 2.4VARIABLE22058{ 2.9VARIABLE22059y 2.7VARIABLE22060@ 2.0VARIABLE22061 3.0VARIABLE22062d- 2.2VARIABLE22063h 2.9VARIABLE22064 2.9VARIABLE22065 3.1VARIABLE22066m 3.0VARIABLE22067pY 2.7VARIABLE22068} 2.2VARIABLE22069 2.5VARIABLE22070t 3.2VARIABLE22071xh 2.8VARIABLE22072 2.5VARIABLE22073 2.8VARIABLE22074 2.9VARIABLE22075 3.0VARIABLE22076 2.8VARIABLE220778 3.0VARIABLE22078B\ 2.9VARIABLE22079FH 2.6VARIABLE22080˵ 2.4VARIABLE22081q 2.4VARIABLE22082K 2.7VARIABLE22083O 2.7VARIABLE22084& 3.0VARIABLE22085 3.4VARIABLE22086S 3.1VARIABLE22087W 2.3VARIABLE22088Z0 3.0VARIABLE22089Q 2.5VARIABLE220901{ 2.6VARIABLE220915y 3.0VARIABLE220929 2.6VARIABLE22093 2.3VARIABLE22094- 2.7VARIABLE22095 3.0VARIABLE22096 2.9VARIABLE22097  2.9VARIABLE22098 2.5VARIABLE22099Y 2.8VARIABLE22100} 3.3VARIABLE23101 2.7VARIABLE23102 3.0VARIABLE23103h 2.9VARIABLE23104 3.0VARIABLE23105 3.0VARIABLE23106 2.5VARIABLE23107 2.9VARIABLE23108 2.5VARIABLE231098 3.6VARIABLE23110\ 3.2VARIABLE23111H 2.7VARIABLE23112 3.0VARIABLE23113.q 2.5VARIABLE23114| 2.8VARIABLE23115 3.2VARIABLE23116& 3.0VARIABLE23117$ 3.8VARIABLE23118^ 2.6VARIABLE23119 2.2VARIABLE231200 3.2VARIABLE23121Q 2.8VARIABLE23122>{ 2.8VARIABLE23123y 2.7VARIABLE23124 3.3VARIABLE23125: 3.2VARIABLE23126=- 2.8VARIABLE23127b 3.0VARIABLE23128f 2.8VARIABLE23129 3.0VARIABLE23130 2.8VARIABLE23131jY 3.8VARIABLE23132n} 2.8VARIABLE23133r 2.8VARIABLE23134 2.6VARIABLE23135~h 3.0VARIABLE23136v 3.4VARIABLE23137z 3.1VARIABLE23138[ 3.0VARIABLE23139 3.1VARIABLE23140 3.1VARIABLE231418 3.1VARIABLE23142\ 2.7VARIABLE23143{H 3.2VARIABLE23144D 3.3VARIABLE23145Hq 3.0VARIABLE23146 2.5VARIABLE23147} 3.0VARIABLE23148M& 3.4VARIABLE23149Q 3.0VARIABLE23150 1.4VARIABLE3/1001\ 1.4VARIABLE3/1002T0 1.3VARIABLE3/1003XQ 1.5VARIABLE3/1004{ 1.4VARIABLE3/1005y 1.7VARIABLE3/10063 1.4VARIABLE3/10077 1.5VARIABLE3/1008- 1.4VARIABLE3/1009 1.5VARIABLE3/1010 1.5VARIABLE3/1011 1.6VARIABLE3/1012 1.4VARIABLE3/1013Y 1.1VARIABLE3/1014} 1.2VARIABLE3/1015 1.5VARIABLE3/1016 1.3VARIABLE3/1017h 1.4VARIABLE3/1018 1.7VARIABLE3/1019 1.5VARIABLE3/1020 1.7VARIABLE3/1021 1.5VARIABLE3/1022 1.0VARIABLE3/10238 1.7VARIABLE3/1024\ 1.9VARIABLE3/1025H 1.6VARIABLE3/1026 1.6VARIABLE3/1027q 1.5VARIABLE3/1028 1.4VARIABLE3/1029( 1.6VARIABLE3/1030& 1.6VARIABLE3/1031 1.5VARIABLE3/1032 1.5VARIABLE3/1033) 1.4VARIABLE3/1034-0 1.5VARIABLE3/1035Q 1.2VARIABLE3/1036{ 1.3VARIABLE3/1037%y 1.4VARIABLE3/1038׌ 1.3VARIABLE3/1039 4 1.5VARIABLE3/1040- 1.3VARIABLE3/1041@ 1.3VARIABLE3/1042 1.3VARIABLE3/1043d 1.6VARIABLE3/1044h 1.9VARIABLE3/1045Y 1.4VARIABLE3/1046} 1.6VARIABLE3/1047l 1.4VARIABLE3/1048p 1.5VARIABLE3/1049h 1.4VARIABLE3/1050 4.7VARIABLE3/2051t 4.5VARIABLE3/2052x 4.9VARIABLE3/2053 4.0VARIABLE3/2054 4.6VARIABLE3/20558 4.5VARIABLE3/2056\ 4.7VARIABLE3/2057H 3.3VARIABLE3/2058ص 4.6VARIABLE3/2059Bq 3.9VARIABLE3/2060F 3.5VARIABLE3/2061˝ 4.2VARIABLE3/2062& 4.0VARIABLE3/2063K 4.7VARIABLE3/2064O 3.6VARIABLE3/2065 4.4VARIABLE3/20660 4.5VARIABLE3/2067Q 4.1VARIABLE3/2068V{ 4.5VARIABLE3/2069Zy 3.9VARIABLE3/2070 4.8VARIABLE3/20711 4.0VARIABLE3/20725- 4.9VARIABLE3/20739 4.7VARIABLE3/2074 4.3VARIABLE3/2075 4.4VARIABLE3/2076 4.8VARIABLE3/2077Y 5.0VARIABLE3/2078 } 4.5VARIABLE3/2079 3.5VARIABLE3/2080 3.8VARIABLE3/2081h 3.7VARIABLE3/2082 3.9VARIABLE3/2083 5.1VARIABLE3/2084 4.5VARIABLE3/2085 4.5VARIABLE3/2086 4.7VARIABLE3/20878 4.4VARIABLE3/2088\ 4.1VARIABLE3/2089H 4.0VARIABLE3/2090 4.4VARIABLE3/2091q 4.6VARIABLE3/2092 4.0VARIABLE3/2093 3.3VARIABLE3/2094.& 4.2VARIABLE3/2095| 4.2VARIABLE3/2096 4.2VARIABLE3/2097 4.3VARIABLE3/2098!0 3.0VARIABLE3/2099;Q 4.1VARIABLE3/2100{ 6.0VARIABLE3/3101y 5.1VARIABLE3/3102 5.9VARIABLE3/3103> 5.6VARIABLE3/3104- 5.8VARIABLE3/3105 6.6VARIABLE3/3106: 4.5VARIABLE3/3107= 6.3VARIABLE3/3108b 5.8VARIABLE3/3109eY 6.1VARIABLE3/3110i} 5.1VARIABLE3/3111 5.3VARIABLE3/3112j 5.5VARIABLE3/3113nh 5.0VARIABLE3/3114r 5.1VARIABLE3/3115 5.3VARIABLE3/3116~ 5.5VARIABLE3/3117v 6.7VARIABLE3/3118z 6.9VARIABLE3/31198 5.0VARIABLE3/3120\ 5.7VARIABLE3/3121H 4.9VARIABLE3/3122 6.7VARIABLE3/3123q 4.9VARIABLE3/3124{ 5.7VARIABLE3/3125D 6.0VARIABLE3/3126H& 4.8VARIABLE3/3127 4.9VARIABLE3/3128} 5.6VARIABLE3/3129M 5.8VARIABLE3/3130P0 6.1VARIABLE3/3131Q 6.4VARIABLE3/3132{ 5.6VARIABLE3/3133Ty 5.1VARIABLE3/3134X 5.6VARIABLE3/3135 6.1VARIABLE3/3136- 5.6VARIABLE3/31373 5.5VARIABLE3/31387 4.8VARIABLE3/3139 5.4VARIABLE3/3140 5.6VARIABLE3/3141 Y 5.1VARIABLE3/3142 } 5.1VARIABLE3/3143  5.9VARIABLE3/3144  5.7VARIABLE3/3145 h 5.2VARIABLE3/3146  5.0VARIABLE3/3147  5.2VARIABLE3/3148  5.4VARIABLE3/3149  5.1VARIABLE3/3150  0.2VARIABLE4a1001 8 0.2VARIABLE4a1002 \ 0.2VARIABLE4a1003 H 0.2VARIABLE4a1004  0.2VARIABLE4a1005 q 0.4VARIABLE4a1006  0.3VARIABLE4a1007 0.2VARIABLE4a1008 & 0.2VARIABLE4a1009  0.1VARIABLE4a1010 ( 0.2VARIABLE4a1011  0.2VARIABLE4a1012 0 0.1VARIABLE4a1013 Q 0.1VARIABLE4a1014 *{ 0.2VARIABLE4a1015 -y 0.4VARIABLE4a1016 0.4VARIABLE4a1017 0.3VARIABLE4a1018 %- 0.3VARIABLE4a1019 ׫ 0.3VARIABLE4a1020  0.2VARIABLE4a1021  0.4VARIABLE4a1022 @ 0.2VARIABLE4a1023 "Y 0.5VARIABLE4a1024 c} 0.2VARIABLE4a1025 g 0.2VARIABLE4a1026 0.4VARIABLE4a1027 h 0.2VARIABLE4a1028 l 0.2VARIABLE4a1029 p 0.2VARIABLE4a1030 0.2VARIABLE4a1031  0.4VARIABLE4a1032 t 0.1VARIABLE4a1033 w8 0.2VARIABLE4a1034 \ 0.2VARIABLE4a1035 H 0.2VARIABLE4a1036 0.2VARIABLE4a1037 q 0.1VARIABLE4a1038 0.2VARIABLE4a1039 ؝ 0.2VARIABLE4a1040 B& 0.3VARIABLE4a1041 F 0.3VARIABLE4a1042  0.2VARIABLE4a1043  0.6VARIABLE4a1044 J0 0.4VARIABLE4a1045 NQ 0.3VARIABLE4a1046  R{ 0.2VARIABLE4a1047 y 0.2VARIABLE4a1048 0.2VARIABLE4a1049 V 0.2VARIABLE4a1050 Z- 1.4VARIABLE4a2051 1.5VARIABLE4a2052 1 1.5VARIABLE4a2053 5 1.3VARIABLE4a2054 9 1.5VARIABLE4a2055 Y 1.3VARIABLE4a2056} 1.6VARIABLE4a2057 1.0VARIABLE4a2058 1.3VARIABLE4a2059 h 1.4VARIABLE4a2060 1.0VARIABLE4a2061 1.5VARIABLE4a2062 1.0VARIABLE4a2063 1.4VARIABLE4a2064 1.3VARIABLE4a20658 1.4VARIABLE4a2066\ 1.5VARIABLE4a2067H 1.0VARIABLE4a2068 1.5VARIABLE4a2069q 1.1VARIABLE4a2070 1.8VARIABLE4a2071 1.3VARIABLE4a2072& 1.5VARIABLE4a2073 1.2VARIABLE4a2074 1.3VARIABLE4a2075. 1.4VARIABLE4a2076+0 1.4VARIABLE4a2077Q 1.7VARIABLE4a2078{ 1.5VARIABLE4a2079!y 1.0VARIABLE4a2080; 1.1VARIABLE4a2081 1.0VARIABLE4a2082- 1.2VARIABLE4a2083 1.6VARIABLE4a2084> 1.5VARIABLE4a2085 1.6VARIABLE4a2086 1.5VARIABLE4a2087`Y 1.3VARIABLE4a2088'} 1.3VARIABLE4a2089a 1.3VARIABLE4a2090e 1.2VARIABLE4a2091ih 1.4VARIABLE4a2092 1.2VARIABLE4a2093j 1.0VARIABLE4a2094n 1.3VARIABLE4a2095r 1.2VARIABLE4a2096 1.3VARIABLE4a20978 1.3VARIABLE4a2098u\ 1.1VARIABLE4a2099yH 1.3VARIABLE4a2100ڵ 2.5VARIABLE4a3101q 1.9VARIABLE4a3102 2.1VARIABLE4a3103 1.8VARIABLE4a3104& 2.2VARIABLE4a3105{ 2.1VARIABLE4a3106D 1.7VARIABLE4a3107H 1.8VARIABLE4a31080 1.8VARIABLE4a3109Q 2.5VARIABLE4a3110L{ 2.0VARIABLE4a3111Py 1.9VARIABLE4a3112 2.1VARIABLE4a3113 2.0VARIABLE4a3114T- 2.4VARIABLE4a3115X 2.3VARIABLE4a3116 1.8VARIABLE4a3117 2.2VARIABLE4a31183 2.3VARIABLE4a31196Y 1.5VARIABLE4a3120} 2.3VARIABLE4a3121 2.0VARIABLE4a3122 2.0VARIABLE4a3123h 1.8VARIABLE4a3124 2.1VARIABLE4a3125 1.8VARIABLE4a3126 1.8VARIABLE4a3127 1.8VARIABLE4a3128 2.1VARIABLE4a31298 1.6VARIABLE4a3130\ 1.9VARIABLE4a3131 H 2.0VARIABLE4a3132 2.2VARIABLE4a3133q 1.5VARIABLE4a3134 1.4VARIABLE4a3135 2.3VARIABLE4a3136& 2.4VARIABLE4a3137 1.8VARIABLE4a3138 1.8VARIABLE4a3139 2.1VARIABLE4a31400 2.4VARIABLE4a3141