C Ohjelma leimauskilpailun tulosten laskentaan C C:\LEIMKIL -alihakemistossa PUUSTO.DAT -tiedosto sek„ C leimattavien puiden numerot k„sitt„v„ tiedosto. C PUUSTO.DAT -tiedostossa koealamittauksen tulokset (formaatti: C metsikk”koealalomake-hyyti„l„ /MARV) C Ohjelma laskee leimaajan kilpailutuloksen, leimatun puuston C puustotunnukset sek„ kasvatettavan puuston puustotunnukset. C Tulostaa PUUSTO.HAR tiedoston - puutiedot HARSI formaatissa C Tulostaa *.LEI tiedoston, jossa leimaus HARSIa varten C Tulostaa *.VIR tiedoston, jossa virheet (j„tt”- ja poistovirheet) C Kaikki tapahtuu C:\LEIMKIL\ hakemistossa C Kutsuttavat aliohjelmat C PITK REGD6 CV2K CV3K RESULT LLASKU CRKOSA (sis„lt„„ funktioita C TTEKO TILOS VALTAP & subroutine aliohj.) INCLUDE 'fgraph.fi' INCLUDE 'fgraph.fd' PARAMETER (MAXP=500) INTEGER X(MAXP),Y(MAXP),JAKSO(MAXP),IPL(MAXP),IPLK(MAXP), +KAYTTO(MAXP),JNRO,NRO(MAXP),IALKU,JAKSOT,LEIMA(MAXP),MALLI(MAXP), +LKML,LKMM,IPLL,MASK,IMIT,ID6MIN,ILASK INTEGER*2 status REAL D(MAXP),DK(MAXP),D6(MAXP),H(MAXP),ID(MAXP),B(MAXP), +IH(MAXP),IKA(MAXP),B10,B11,B20,B21,PIT(MAXP),TIL(MAXP),A,BE,PALA, +VP,VJ,TULOS,PPPA,PDGM,PHGM,PTIL,JPPA,JDGM,JTIL,JHGM,EPPA,EDGM, +ETIL,EHGM,VTAB(5),VTUK(MAXP),VKUI(MAXP),VHUK(MAXP),VTOT(MAXP), +HPTAB(3),SUMTUK,SUMKUI,SUMHUK,SUMTOT,SUMTIL,A1,RLUKU,RTUK,RKUI, +RPOI,SUMD2,SUMD3,DTUK,D6TUK,HTUK,TUKKES,VPLAJI(4),VOSUUS(4), +HDOM,JHDOM CHARACTER CMALLI*12,APU1*23,LUKU*32,RIVI1*30,NIMI*30,CLEIM*5, +CLEIMA*12,APU2*23,CMAL*6,VAST*1,CNUM*3 RECORD /rccoord/ curpos OPEN (15,FILE='C:\LEIMAUS\MALLI.DAT',STATUS='OLD') OPEN (16,FILE='C:\LEIMAUS\PUUSTO.DAT',STATUS='OLD') status = setvideomode ($MAXRESMODE) CALL clearscreen ($GCLEARSCREEN) CALL settextposition (3,1,curpos) status = settextcolor(2) CALL outtext ('Enter FILENAME for the trees to be removed:') READ (5,'(A12)') CLEIMA APU2(1:23)='C:\LEIMAUS\'//CLEIMA(1:12) OPEN (14,FILE=APU2,STATUS='OLD') C puusto.dat -tiedoston tietuerakenne (kanava 16): C 1- 3 puun numero NRO Integer C 4- 6 x-koordinaatti,dm IX Integer C 7- 9 y-koordinaatti,dm IY Integer C 10 jakso 1 tai 2 JAKSO Integer C 11 puulaji 1,2,..,9 IPL Integer C 12 puuluokka IPLK Integer C 13-15 d1.3 D Real C 16 k„ytt” KAYTTO Integer C 17-18 dk DK Real C 19-20 d6 D6 Real C 21-23 h H Real C 24-25 id ID Real C 26-27 2xB B Real C 28-29 ih IH Real C 30-32 ik„ IKA Real C *************************************************** C puutietojen luku C:\LEIMKIL\PUUSTO.DAT -tiedostosta CALL settextposition (5,1,curpos) DO 3 I=1,192 READ (16,'(A32)',END=5) LUKU IF (LUKU(13:13).EQ.'0') GOTO 3 READ(LUKU,201)NRO(I),X(I),Y(I),JAKSO(I),IPL(I),IPLK(I),D(I), +KAYTTO(I),DK(I),D6(I),H(I),ID(I),B(I),IH(I),IKA(I) 201 FORMAT(I3,I3,I3,I1,I1,I1,F3.0,I1,F2.0,F2.0,F3.0,F2.0,F2.0,F2.0, +F3.0) D(I)=D(I)/10. H(I)=H(I)/10. WRITE(CNUM,'(I3)') I CALL outtext (cnum) CALL settextposition (5,1,curpos) 3 CONTINUE 5 CONTINUE I=192 C I on kokonaislukulaskuri puiden lukum„„r„lle (192) C ***************************************************************** C pituusk„yr„n parametrien B10 ja B11 laskenta 1. jakson puille C h=B10(d/(d+B11)) aliohjelmalla PITK C ! ongelma ! jos ei muita kuin koepuita J1NRO saa v„„r„n arvon J1NRO=192 C WRITE (6,*) ' The number of sample trees is ', 192 CALL settextposition (8,1,curpos) CALL outtext ('Computing height curves, please waite!') IALKU=1 CALL PITU(IALKU,J1NRO,D,H,B10,B11) c B10=47.00 c B11=32.10 c WRITE (6,*) ' Parameters have been calculated ' c WRITE (6,*) ' 1 layer B0 ',B10,' B1 ',B11 C (Etsint„algoritmi perustuu) PUUSTO.DAT tiedoston tietueiden j„rjestykseen: C 1. jakson koepuut, 2. jakson koepuut, 1. jakson lukupuut, 2. jakson lukupuut C Lasketaan pituudet kaikille puille taulukkoon PIT DO 12 M=1,I PIT(M)=1.3 + (D(M)**2) / (B10 + B11*D(M)) **2 12 CONTINUE C *********************************************************** C Lasketaan lineaarinen regressio d1.3 ja d6 v„lille aliohjelmalla REGD6 C Tarkistetaan ensin, ett„ D6 on mitattu 'riitt„v„n monesta koepuusta' C 'Riitt„v„n monta' ID6MIN=6 CALL REGD6(D,D6,I,A,BE) C D6 yleist„minen lukupuille ohitetaan, jollei minimim„„r„ t„yty C Annetaan d6 kaikille lukupuille, koepuiden d6 ennallaan DO 17 M=1,I IF(D6(M).LT.0.1 .AND. PIT(M).GT.7.6) THEN D6(M)=A+(BE*D(M)) END IF 17 CONTINUE C ****************************************************************** C Lasketaan tilavuudet kaikille koealan puille, regressiomalleilla C aliohjelmat CV2K ja CV3K DO 13 M=1,I IF (D6(M).LT.0.1 .OR. PIT(M).LT.7.6) GOTO 14 TIL(M)=CV3K(IPL(M),D(M),D6(M),PIT(M)) 14 CONTINUE TIL(M)=CV2K(IPL(M),D(M),PIT(M)) 13 CONTINUE C ***************************************************************** C Luetaan LEIMAAJAN poistettavaksi tarkoittamien puiden numerot C kanavalta 14 (tiedosto avattuna) 20 CONTINUE READ (14,'(A30)') RIVI1 IF (RIVI1(1:1).EQ.'1' .OR. RIVI1(1:1).EQ.'2' .OR. + RIVI1(1:1).EQ.'3' .OR. RIVI1(1:1).EQ.'4' .OR. + RIVI1(1:1).EQ.'5' .OR. RIVI1(1:1).EQ.'6' .OR. + RIVI1(1:1).EQ.'7' .OR. RIVI1(1:1).EQ.'8' .OR. + RIVI1(1:1).EQ.'9' .OR. RIVI1(1:1).EQ.'0' .OR. + RIVI1(1:1).EQ.' ') THEN CALL clearscreen ($GCLEARSCREEN) CALL settextposition (4,1,curpos) CALL outtext (' file ') CALL outtext (APU1) CALL outtext (' 1. line should read your NAME') GOTO 998 END IF NIMI(1:30)=RIVI1(1:30) 21 CONTINUE LKML=LKML+1 READ (14,'(A5)',ERR=998) CLEIM IF (CLEIM(1:5).EQ.'LOPPU') GOTO 22 READ (CLEIM,'(I5)',ERR=999) LEIMA(LKML) IF(LEIMA(LKML).GT.999) GOTO 998 GOTO 21 22 CONTINUE LKML=LKML-1 CALL clearscreen ($GCLEARSCREEN) status = rectangle ($GBORDER,100,100,400,400) DO 9091 ib=1,192 newx = x (ib) + 100 newy = - ((y(ib)-200) - 200) newx1 = newx - INT (d(ib))/2 newx2 = newx + INT (d(ib))/2 newy1 = newy - INT (d(ib))/2 newy2 = newy + INT (d(ib))/2 DO 9092 ic=1,LKML IF (ib.EQ.LEIMA(ic)) THEN status = ellipse($GFILLINTERIOR,newx1,newy1,newx2,newy2) GOTO 9091 END IF 9092 CONTINUE status = ellipse ($GBORDER,newx1,newy1,newx2,newy2) 9091 CONTINUE READ (*,*) CALL clearscreen ($GCLEARSCREEN) C LKML=leimattujen puiden lukum„„r„ C *********************************************************** C Luetaan MALLILEIMAUS ja C:\LEIMAUS\ -hakemistosta C kanavalta 15, PALA=0.09 23 CONTINUE PALA=0.09 24 CONTINUE LKMM=LKMM+1 READ (15,*,END=26) MALLI(LKMM) IF (MALLI(LKMM).GT.999 .OR. MALLI(LKMM).EQ.0) GOTO 997 GOTO 24 26 CONTINUE LKMM=LKMM-1 C LKMM = mallileimattujen puiden lukum„„r„ C ****************************************************************** C T„ss„ vaiheessa on taulukossa TIL puiden tilavuudet I kpl C LEIMA leimaajan poistamien numerot LKML(kpl) C MALLI mallileimauksessa poistetut LKMM(kpl) C D puiden l„pimitat, cm I kpl C PIT puiden pituudet, m I kpl C D6 yl„l„pimitat, cm I kpl C TIL tilavuudet, dm3 CALL RESULT(LEIMA,LKML,MALLI,LKMM,I,VP,VJ,TULOS,APU2) C RESULT laskee tuloksen ja kirjoitta c:\leimaus\*.vir C -tiedoston, jossa alkuosa * on sama kuin leimaustiedostossa C Tiedosto sis„lt„„ virhepuiden numerot, sek„ J tai P kirjaimen C J v„„rin j„tetyille; P v„„rin poistetuille C VP v„„rin poistettujen lukum„„r„ (poistovirhe) C VJ v„„rin j„tettyjen lukum„„r„ (j„tt”virhe) C TULOS leimaajan tulos (puiden lukum„„r„ v„hennettyn„ edelliset virheet) CALL LLASKU(LEIMA,LKML,MALLI,LKMM,D,D6,PIT,TIL,I,PPPA,PTIL,PDGM, +PHGM,JPPA,JTIL,JDGM,JHGM,EPPA,ETIL,EDGM,EHGM,PALA) CALL VALTAP(D,I,PALA,B10,B11,HDOM,JHDOM,LEIMA,LKML) C PPPA poistuman pohjapinta-ala, PTIL poistuman tilavuus, C PDGM poistuman keskil„pimitta, PHGM poistuman keskipituus, C JPPA j„„v„n puuston pohjapinta-ala, JTIL j„„v„n puuston tilavuus, C JDGM j„„v„n puuston keskil„pimitta, JHGM j„„v„n puuston keskipituus, C EPPA pohjapinta-ala ennen harvennusta, ETIL tilavuus ennen harvennusta, C EDGM keskil„pimitta ennen harvennusta, EHGM keskipituus ennen harvennusta C PALA Koealan pinta-ala (ha) C HDOM Valtapituus JHDOM J„„v„n puuston valtapit. CALL clearscreen ($GCLEARSCREEN) CALL settextposition (4,1,curpos) CALL outtext( ' Do you wish to run HARSI-program later y/n:') READ (5,'(A1)') VAST IF (VAST(1:1).EQ.'N' .OR. VAST(1:1).EQ.'n') GOTO 27 CALL TTEKO(NRO,X,Y,IPL,D,PIT,I,APU2,LEIMA,LKML) C TTEKO aliohjelma tulostaa *.har ja *.lei tiedostot harsia varten C C:\HARSI\ -alihakemistoon 27 CONTINUE C Lasketaan VTAB(1:5) REAL-taulukkoon C VTAB(1) Tukkiosan tilavuus [dm3] C VTAB(2) Kuituosan tilavuus [dm3] C VTAB(3) Latvahukkaosan tilavuus [dm3] C VTAB(4) Mahd. tyveys (ei koske t„t„ ohjelmaa) C VTAB(5) Rungon tilavuus [dm3] IMIT=1 MASK=0 DO 28 K=1,I IPLL=IPL(K) IF (IPL(K).EQ.4) IPLL=3 IF (IPL(K).EQ.5) IPLL=3 IF (IPL(K).EQ.6) IPLL=3 IF (IPL(K).EQ.7) IPLL=3 IF (IPL(K).EQ.8) IPLL=2 CALL CRKOSA(IPLL,IMIT,D(K),D6(K),PIT(K),MASK,HPTAB,VTAB) VTUK(K)=VTAB(1) VKUI(K)=VTAB(2) IF(IPLK(K).EQ.2) VKUI(K)=VKUI(K)+VTUK(K) VHUK(K)=VTAB(3) VTOT(K)=VTAB(5) 28 CONTINUE C Lasketaan koealakohtaiset summatunnukset, tukki, kuitu, hukka ja kokonais C tilavuudelle, runkok. ja regressiomallilla DO 29 K=1,I SUMTUK=SUMTUK+VTUK(K) SUMKUI=SUMKUI+VKUI(K) SUMHUK=SUMHUK+VHUK(K) SUMTOT=SUMTOT+VTOT(K) SUMTIL=SUMTIL+TIL(K) 29 CONTINUE C Hehtaarikohtaisiksi A1=1./PALA SUMTUK=SUMTUK*A1 SUMHUK=SUMHUK*A1 SUMKUI=SUMKUI*A1 SUMTOT=SUMTOT*A1 SUMTIL=SUMTIL*A1 RLUKU=REAL(I)*A1 DO 30 K=1,I IF(IPLK(K).EQ.1) RTUK=RTUK+1. IF(IPLK(K).EQ.2) RKUI=RKUI+1. IF(IPLK(K).EQ.4) RPOI=RPOI+1. 30 CONTINUE RTUK=RTUK*A1 RKUI=RKUI*A1 RPOI=RPOI*A1 C Tukkiosan keskitilavuus TUKKES, haetaan keskim. tukkirungon l„pimitta, C sille d6 ja pituus ja kutsutaan CRKOSAA DO 31 K=1,I IF(IPLK(K).EQ.1) THEN SUMD2=SUMD2+D(K)**2 SUMD3=SUMD3+D(K)**3 END IF 31 CONTINUE DTUK=SUMD3/SUMD2 D6TUK=A+(BE*DTUK) HTUK=1.3 + (DTUK**2)/(B10+B11*DTUK**2) WRITE (6,*) ' Dominant tree species ? (1,2,3,..,9)>' READ (5,'(I1)')IIPL CALL CRKOSA (IIPL,IMIT,DTUK,D6TUK,HTUK,MASK,HPTAB,VTAB) TUKKES=VTAB(1) C Lasketaan tilavuudet puulajeittain ja lasketaan n„ist„ puulajiosuudet C Laskentaryhmin„ m„nty,kuusi,koivut,muulp CALL TILOS (IPL,VTOT,I,VPLAJI,VOSUUS,PALA) CALL TEXT (VPLAJI,VOSUUS,SUMTUK,SUMKUI,SUMHUK,SUMTOT,SUMTIL, +RLUKU,RTUK,RKUI,RPOI,TUKKES,PPPA,PDGM,PTIL,PHGM,JPPA,JDGM, +JTIL,JHGM,EPPA,EDGM,ETIL,EHGM,PALA,D,D6,PIT,TIL,VTUK, +VKUI,VHUK,VTOT,I,TULOS,VP,VJ,LKML,HDOM,JHDOM) GOTO 999 997 WRITE (6,*) ' There is something wrong with ',APU2 GOTO 999 998 WRITE (6,*) ' There is something wrong with ',APU1 999 CONTINUE END C ****** C AITKEN C ****** SUBROUTINE AITKEN(F,PO,X,YY,M,Z) DIMENSION F(M),X(M),PO(M,M) DO 100 K=1,M PO(K,1)=F(K) 100 CONTINUE DO 200 I=1,M-1 DO 200 K=I+1,M PO(K,I+1)=((X(K)-YY)*PO(I,I)-(X(I)-YY)*PO(K,I))/(X(K)-X(I)) 200 CONTINUE Z=PO(M,M) RETURN END C *************** C AKANT1 FUNCTION C *************** FUNCTION AKANT1(IP,D,H) C---KANNON PITUUS (M) C---MAT.OS/JTH 16.11.77 C---PARAMETRIARVOJEN TESTAUS LISATTY 30.6.79 (JTH) C---LEHTIKUUSI LIS[TTY 7.10.1982 / C-G S C-- IP = PUULAJI, 1<=IP<=7 C-- D = KUORELLINEN RINNANKORK. LAPIMITTA (CM), 0.1<=D<=85. C-- H = PUUN PITUUS (M), 1<=H<=40. C-- FUNKTION ARVO >= 0.07, PAITSI JOS JOKIN YO EHDOISTA EI TOTEUDU, C JOLLOIN ARVO = 0. IF(IP.LT.1.OR.IP.GT.7) GO TO 5 IF(D.LT.0.1.OR.D.GT.85.) GO TO 5 IF(H.LT.1.0.OR.H.GT.40.) GO TO 5 GO TO (1,2,3,3,8,7,6) IP C---LEHTIKUUSI 6 G=0.803156*H+0.101395*D GO TO 4 C---LEPPA 7 G=.7168*D+.016653*H GO TO 4 C---HAAPA 8 GA=.509136+.840372*ALOG(D) G=EXP(GA)-1. GO TO 4 C---MANTY 1 G=.09522*H+.4456*D GO TO 4 C---KUUSI 2 G=.56*H+.5089*D GO TO 4 C---KOIVU 3 G=.497936*H+.4862*D 4 AKANT1=G/100. RETURN 5 AKANT1=0. WRITE(5,100) 100 FORMAT(' ALIOHJELMAN AKANT1 ARGUMENTEISSA VIRHE, SAADUT + ARVOT:') WRITE(5,101)IP,D,H 101 FORMAT(' PUULAJI =',I2,' D =',F5.1,' H =',F5.1) RETURN END C ********** C * CPOLY3 * C ********** C C OHJELMA LASKEE KORJAUSPOLYNOMIN Y = B(1)*X + B(2)*X**2 + B(3)*X**3 C KERTOIMET SITEN, ETT[ POLYNOMI KULKEE SEURAAVAN NELJ[N PISTEEN C KAUTTA: C X1 = P(1) = (H-1.3)/H (RINNANKORKEUS) Y(X1) = 0 C X2 = P(2) = (H-6)/H (6 M) Y(X2) = P(3) = DIF C X3 = P(4) (TUKIPISTE) Y(X3) = P(5) C X4 = 0 (LATVA) Y(X4) = 0 C MATEMAATTINEN OSASTO - RUNKOK[YR[T 28.8.1980 / C-G S C C SUBROUTINE CPOLY3 (P,B) C DIMENSION P(5), B(3) C CON1 = P(3) / (P(2) * (P(2)-P(1))) CON2 = P(5) / (P(4) * (P(4)-P(1))) C B(3) = (CON1-CON2) / (P(2)-P(4)) B(2) = CON1 - B(3) * (P(1)+P(2)) B(1) = P(1) * (P(2)*B(3) - CON1) C RETURN END C *** C CRK C *** SUBROUTINE CRK(IPL,D,D6,H,HMIT,C) C---TEKEE RUNKOKAYRAN JA KORJAA SITA ENNUSTEYHTALOIDEN AVULLA C---MAT.OS/JTH 24.4.80 KORJATTU 15.9.1980 / C-G S C IPL = PUULAJI, 1<=IPL<=9 C D = D13 (CM) C D6 = LAPIMITTA (CM) C H = PUUN PITUUS (M) C HMIT= MITTAUSTASON KORKEUS MAANPINNASTA (M),(ESIM C JUURENNNISKA) C C = RUNKOKAYRAN KERTOIMET, D20 % JA H MAANPINNASTA DIMENSION C(10),P(5),B(3) C C---VIRHETILANTEET IVIR=0 IF(IPL.LT.1.OR.IPL.GT.9)IVIR=1 IF(D.LT.0.1.OR.D.GT.90.)IVIR=1 IF(D6.LT.0..OR.D6.GT.D+1.)IVIR=1 IF(H.LT.1.4.OR.H.GT.40.)IVIR=1 IF(HMIT.LT.0..OR.HMIT.GT.1.)IVIR=1 IF(IVIR.GT.0) GO TO 90 C C---APUSUUREET, PERUSKERTOIMET , C(9) JA C(10) DA=D H13=1.3 IF(HMIT.GE.0.01.AND.H.GE.1.5)H13=H13+HMIT CALL CRKPK(IPL,C) C(10)=H+HMIT C(9)=D/CRKT(H13,C) IF(D6.GT.0..AND.H.GT.6.1) GO TO 1 C C---KORJAUSVEKTORI P D:N JA H:N PERUSTEELLA KUN IPL<=4. IF(IPL.GT.4) GO TO 91 IF(HMIT.GE.0.01) DA=C(9)*CRKT(1.3,C) CALL CRKKD(IPL,DA,C,P) GO TO 9 C C---KORJAUSVEKTORI P D:N, D6:N JA H:N PERUSTEELLA, IPL<=9 1 H6=6.+HMIT D6A=D6 IF(HMIT.LT.0.01) GO TO 2 DA=C(9)*CRKT(1.3,C) D6A=D6/CRKT(H6,C)*CRKT(6.,C) 2 P(1)=1.-H13/C(10) P(2)=1.-H6/C(10) P(3)=D6/C(9)-CRKT(H6,C) CALL CRKKD6(IPL,DA,D6A,HMIT,C,P) C C---RUNKOKAYRAN KERTOIMIEN MUUTOKSET 9 CALL CPOLY3(P,B) C DO 8 I=1,3 8 C(I)=C(I)+B(I) C(9)=D/CRKT(H13,C) GO TO 91 C C---VIRHEILMOITUS 90 WRITE(5,60)IPL,D,D6,H,HMIT 60 FORMAT(' ALIOHJELMAN RK PARAMETREISSA IVIRHE, SAADUT ARVOT:'/ +' IPL=',I3,' D =',G5.1,' D6 =',G5.1,' H =',G5.1,' HMIT=',G5.1) DO 89 I=1,9 89 C(I)=0. C C---PALUU OHJELMASTA 91 RETURN END C ******** C * CRKD * C ******** C C FUNKTIO LASKEE RUNKOK[YR[LT[ PUUN L[PIMITAN (CM) KORKEUDELLA HX C (M) C MAANPINNAN TASOSTA. C(1)-C(8) OVAT RUNKOK[YR[N KERTOIMET, C(9) ON C L[PIMITTA (CM) KORKEUDELLA 0.2*H JA C(10) ON PUUN PITUUS H (M). C MATEMAATTINEN OSASTO - RUNKOK[YR[T 1.9.1980 / C-G S C C REAL FUNCTION CRKD (HX,C) C DIMENSION C(10) C IF (HX .LT. 0.) GO TO 1 IF (HX .GT. C(10)) GO TO 2 C C FUNKTIO CRKT LASKEE SUHTEEN CRKD/C(9): C CRKD = C(9) * CRKT(HX,C) RETURN C 1 WRITE (5,601) HX 601 FORMAT (' CRKD: KORKEUS',F8.2,' M < 0') GO TO 3 2 WRITE (5,602) HX, C(10) 602 FORMAT (' CRKD: KORKEUS',F8.2,' M > PUUN PITUUS',F7.2,' M') 3 CRKD = 0. RETURN END C **** C CRKH C **** FUNCTION CRKH(DX,C) C---ETSII L[PIMITTAA DX VASTAAVAN KORKEUDEN RK:LT[ C---MAT.OS/JTH 20.4.80 8.10.1980 / C-G S DIMENSION C(10) H=C(10) G=1. HYPPY=1. IF(0.0.GT.DX.OR.DX.GT.85..OR.1.4.GT.H.OR.H.GT.40.)GO TO 4 L=H/HYPPY HR=H IF(L*HYPPY.EQ.H) HR=H+.01 C---1.HAARUKOINTI DO 1 I=1,L HR=HR-HYPPY IF(C(9)*CRKT(HR,C).GT.DX) GOTO 2 IF(I.NE.L) GO TO 1 G=-1. HR=HYPPY 1 CONTINUE C---TARKENNUS 2 DO 3 I=1,7 HYPPY=HYPPY/2. HR=HR+G*HYPPY DR=C(9)*CRKT(HR,C) IF(DR.GT.DX)G=1. IF(DR.LE.DX)G=-1. 3 CONTINUE CRKH=HR GO TO 7 C---VIRHEILMOITUKSET 4 WRITE(5,10) 10 FORMAT(' ALIOHJELMAN CRKH ARGUMENTEISSA VIRHE') IF(0.1.LE.DX.AND.DX.LE.85.0)GO TO 5 WRITE(5,11)DX 11 FORMAT(' L[PIM. DX=',F9.2,' OLTAVA V[LILL[ (0.1,85)') 5 IF(1.4.LE.H.AND.H.LE.40.) GO TO 6 WRITE(5,12)H 12 FORMAT(' PITUUS H=',F9.2,' OLTAVA V[LILL[ (1.4,40)') 6 CRKH=0. 7 RETURN END C ********* C * CRKKD * C ********* C C ALIOHJELMA MUODOSTAA P-VEKTORIN D:N JA H:N AVULLA C LASKETTAVALLE C RUNKOK[YR[LLE. KORJAUSPOLYNOMI B ON LASKETTAVA SITEN, ETT[ C B(0) = 0 LATVA C B(P(1)) = 0 0.1*H C B(P(2)) = P(3) 0.4*H C B(P(4)) = P(5) 0.7*H C IPL ON PUULAJI, D ON RINNANKORKEUSL[PIMITTA (CM). C(1)-C(8) OVAT C RUNKOK[YR[N PERUSKERTOIMET JA C(10) PUUN PITUUS H (M). KOSKA B:N C ARVO RINNANKORKEUSL[PIMITTAPISTEESS[ EI OLE 0 VOIDAAN C(9) (D80) C LASKEA VASTA RUNKOK[YR[KORJAUKSEN J[LKEEN! C MATEMAATTINEN OSASTO - RUNKOK[YR[T 23.9.1980 / C-G S C C SUBROUTINE CRKKD (IPL,D,C,P) C DIMENSION C(10) DIMENSION P(5) C C LASKETAAN ALUKSI TARVITTAVAT MUUTTUJAT: C H = C(10) DH = D / (H-1.3) DH2 = DH**2 DL = ALOG(D) HL = ALOG(H) D2 = D**2 C C LASKETAAN PERUSK[YR[N ARVOT T1, T4 JA T7 SEK[ KORJAUKSET Y1, Y4 JA C Y7 C KORKEUKSILLA 0.1*H, 0.4*H JA 0.7*H: C GO TO (1,2,3), IPL C C M[NTY: C 1 T1 = 1.100553 T4 = 0.8585458 T7 = 0.5442665 C Y1 = 0.26222 - 0.0016245*D + 0.010074*H + 0.06273*DH - $ 0.011971*DH2 - 0.15496*HL - 0.45333/H Y4 = -0.38383 - 0.0055445*H - 0.014121*DL + 0.17496*HL + 0.62221/H Y7 = -0.179 + 0.037116*DH - 0.12667*DL + 0.18974*HL GO TO 4 C C KUUSI: C 2 T1 = 1.0814409 T4 = 0.8409653 T7 = 0.4999158 C Y1 = -0.003133*D + 0.01172*H + 0.48952*DH - 0.078688*DH2 - $ 0.31296*DL + 0.13242*HL - 1.2967/H Y4 = -0.0065534*D + 0.011587*H - 0.054213*DH + 0.011557*DH2 + $ 0.12598/H Y7 = 0.084893 - 0.0064871*D + 0.012711*H - 0.10287*DH + $ 0.026841*DH2 - 0.01932*DL GO TO 4 C C KOIVU: C 3 T1 = 1.084544 T4 = 0.8417135 T7 = 0.4577622 C Y1 = 0.59848 + 0.011356*D - 0.49612*DL + 0.46137*HL - $ 0.92116/DH + 0.25182/DH2 - 0.00019947*D2 Y4 = -0.96443 + 0.011401*D + 0.13870*DL + 1.5003/H + $ 0.57278/DH - 0.18735/DH2 - 0.00026*D2 Y7 = -2.1147 + 0.79368*DL - 0.51810*HL + 2.9061/H + $ 1.6811/DH - 0.40778/DH2 - 0.00011148*D2 C 4 Y1 = SIGN(AMIN1(ABS(Y1),0.1),Y1) Y4 = SIGN(AMIN1(ABS(Y4),0.1),Y4) Y7 = SIGN(AMIN1(ABS(Y7),0.1),Y7) C C LASKETAAN P-VEKTORI. KORJAUKSET ON LASKETTU OLETUKSELLA B(P(1)) C = Y1. C KOSKA KUITENKIN B(P(1)):N ON OLTAVA 0 ON KORJAUKSET TEHT[V[ C K[YR[LLE C JOKA SAADAAN KERTOMALLA OLETUSK[YR[[ LUVULLA (T1+Y1)/T1: C P(1) = 0.9 P(2) = 0.6 P(3) = T1/(T1+Y1) * (T4+Y4) - T4 P(4) = 0.3 P(5) = T1/(T1+Y1) * (T7+Y7) - T7 C RETURN END C ********** C * CRKKD6 * C ********** C C ALIOHJELMA M[[R[[ RUNKOK[YR[N KORJAUSPOLYNOMIN APUPISTEEN P(4) C SEK[ C KORJAUSPOLYNOMIN ARVON P(5) APUPISTEESS[ D:N, D6:N JA H:N C PERUSTEELLA. C IPL = PUULAJI C D = RINNANKORKEUSL[PIMITTA (CM) C D6 = KUUDEN METRIN L[PIMITTA (CM) C HMIT = MITTAUSTASON KORKEUS MAANTASALTA (M) C C(1)-C(8) OVAT RUNKOK[YR[N PERUSKERTOIMET, C(9) L[PIMITTA (CM) KOR- C KEUDELLA 0.2*H (D80) JA C(10) PUUN PITUUS H. C MATEMAATTINEN OSASTO - RUNKOK[YR[T 23.9.1980/3.4.1981 / C-G S C PITUUSRAJAT MUUTETTU 10.12.1981 / C-G S C PULLISTUMAKORJAUS 8.1.1982 / C-G S C LEHTIKUUSI LIS[TTY 5.11.1982 / C-G S C C SUBROUTINE CRKKD6 (IPL,D,D6,HMIT,C,P) C DIMENSION C(10) DIMENSION P(5) C DIMENSION F(4), FX(4), PO(4,4) C C LASKETAAN TARVITTAVAT APUSUUREET H, HSD, HD, D6D, AT JA DIF: C H = C(10) HSD = (H-6.) / D6 HD = (H-1.3) / D D6D = (D-D6) / D AT = D6 / C(9) DIF = P(3) IF (HMIT .GE. 0.01) DIF = AT - CRKT(6.,C) C C P(4) JA P(5); P(5):LLE ON JOKAISELLE PUULAJILLE KOLME KAAVAA: C - KAAVAT 11-13 ISOILLE PUILLE > 13.05 M C - KAAVAT 21-23 PIENILLE PUILLE < 8.35 M C - KUN 8.35 M < PUUN PITUUS < 13.05 M LASKETAAN TULOS MOLEMPIEN C EDELL[MAINITTUJEN KAAVOJEN PERUSTEELLA. C IF (H .LT. 8.35) GO TO 20 C C ISOT PUUT: C P(4) = 0.5 - 3./H GO TO (11,12,13,13,13,13,17,11,13), IPL C C M[NTY: C 11 CX = ABS(H-10.7) / 9.4 XX = DIF * (1.+CX) * (H+4.) / (H-0.7) P(5) = -0.31017 + 1.2036*DIF + 0.096447*HSD + 0.32458*D6D - $ 0.16066*XX + 0.21072*AT GO TO 15 C C KUUSI: C 12 P(5) = -0.51243 + 1.0204*DIF - 0.10077*HD + 0.19315*HSD + $ 0.52327*D6D + 0.40615*AT GO TO 15 C C KOIVU: C 13 P(5) = -0.86846 + 0.61482*DIF + 0.079054*HD + 0.84701*D6D + $ 0.71465*AT GO TO 15 C C LEHTIKUUSI: C 17 P(5) = -0.719394 + 1.00321*DIF + 0.0825812*HSD + 0.778782*D6D + $ 0.547997*AT C 15 IF (H .GE. 13.05) RETURN P4ISO = P(4) P5ISO = P(5) C C PIENET PUUT: C 20 P(4) = 1. - 3.65/H GO TO (21,22,23,23,23,23,27,21,23), IPL C C M[NTY: C 21 P(5) = -0.039003 + 0.39739*DIF + 0.042010*HD - 0.07917*HSD + $ 0.062884*AT GO TO 25 C C KUUSI: C 22 P(5) = -0.037635 + 0.53502*DIF + 0.072073*HD -0.065727*HSD + $ 0.032871*AT GO TO 25 C C KOIVU: C 23 P(5) = -0.083783 + 0.37229*DIF + 0.097625*AT GO TO 25 C C LEHTIKUUSI: C 27 CX = ABS(10.7-H) / 9.4 XR = DIF * (1.+CX) XX = 4.7 * XR / (H-1.3) P(5) = 1.10708 - 1.24663*D6D - 0.739426*XX - 0.978947*AT + $ 0.613124*XR C 25 IF (H .LT. 8.35) GO TO 30 C C KESKIKOKOISET PUUT: C FX(1) = P(1) F(1) = 0. FX(2) = P(2) F(2) = DIF FX(3) = P4ISO F(3) = P5ISO FX(4) = 0. F(4) = 0. TT = (H-8.35) / 4.7 IF (H .GT. 10.7) GO TO 26 CALL AITKEN(F,PO,FX,P(4),4,YY) P(5) = TT*YY + (1.-TT)*P(5) GO TO 30 C 26 FX(3) = P(4) F(3) = P(5) P(4) = P4ISO CALL AITKEN(F,PO,FX,P(4),4,YY) P(5) = TT*P5ISO + (1.-TT)*YY C C KORJAUS MAHDOLLISTA PULLISTUMAA VASTAAN V[LILL[ (1.3 M, 6 M): C 30 IF (ABS(DIF) .LT. 0.15) RETURN T13 = D / C(9) XE = CRKT(3.65,C) IF (H .GT. 10.7) GO TO 31 XX = XE + P(5) IF (XX .LT. T13) RETURN P(5) = 0.5*(T13+AT) - XE RETURN C C KORJAUS MAHDOLLISTA PULLISTUMAA VASTAAN KUN 10.7 M < H < 13.05 M: C 31 FX(3) = P4ISO F(3) = P(5) T2T = (H-10.7) / 2.35 XT = 1. - 3.65/H CALL AITKEN(F,PO,FX,XT,4,YY) XX = XE + YY IF (XX .LT. T13) RETURN TX = 0.5*(T13+AT) - XE FX(3) = XT F(3) = TX CALL AITKEN(F,PO,FX,P(4),4,YY) P(5) = T2T*P(5) + (1.-T2T)*YY C RETURN END C CRKOSA C ====== C C Aliohjelma laskee puun tukki-, kuitu- ja latvahukkaosien tilavuudet, C mahdollisen tyveyksen tilavuuden sek{ n{iden osien yhteistilavuuden C puulajin, D:n D6:n ja H:n perusteella. C Osien alku- ja p{{ttymiskorkeudet voidaan antaa ohjelmassa mutta C ellei kaikkia anneta lasketaan puuttuvat seuraavasti: C C Tukkiosan alkukorkeus: estimoitu juurenniska, v{hint{{n 10 cm. C Tukkiosan p{{ttymiskorkeus (= kuituosan alkukorkeus): PMP-systeemin C p|lkytysalgoritmi, C Kuituosan p{{ttymiskorkeus (= latvahukkaosan alkukorkeus): korkeus C jolla (kuoreton) l{pimitta on 6.3 cm (m{nty), 6.5 cm (kuusi, C koivu). C C Sy|tt|parametrit: C C IPL Puulaji 1 = m{nty 2 = kuusi 3 = koivu C IMIT = 1 jos kaikki mitat (D, D6, H, HPTAB) ovat maanpinnasta, C = 2 jos kaikki mitat ovat juurenniskalta. C D Puun rinnankorkeusl{pimitta, cm C D6 Puun kuuden metrin l{pimitta (jos olemassa ja tunnettu), C cm; muussa tapauksessa 0.0 C H Puun pituus, m, maanpinnasta tai juurenniskalta IMIT:n C mukaan C MASK Kolminumeroinen kokonaisluku, abc, joka kertoo, onko C HPTAB:ssa annettu alku/p{{ttymiskorkeuksia vai ei: C a = 1 jos tukkiosan alkukorkeus on annettu, C = 0 ellei ole; C b = 1 jos tukkiosan p{{ttymiskorkeus on annettu, C = 0 ellei ole; C c = 1 jos kuituosan p{{ttymiskorkeus on annettu, C = 0 ellei ole. C Esim. C MASK = 0 Ohjelma laskee kaikki korkeudet C MASK = 110 Ohjelma laskee vain kuituosan p{{ttymiskor- C keuden; tukkiosan alku- ja p{{ttymiskorkeu- C det annetaan HPTAB:ssa. C C Sy|tt|/tulosparametrit: C C HPTAB(1:3) HPTAB(1) Tukkiosan alkukorkeus, m C HPTAB(2) Tukkiosan p{{ttymiskorkeus, m C HPTAB(3) Kuituosan p{{ttymiskorkeus, m C Tulotilanteessa on HPTAB:ssa niin monta korkeutta kuin C MASK-parametrissa on ykk|si{. Paluutilanteessa sis{lt{{ C HPTAB kaikki alku/p{{ttymiskorkeudet (annetut ja lasketut). C C Huom! Jollei ole tukkia on HPTAB(1) = HPTAB(2). Vastaa- C vasti HPTAB(2) = HPTAB(3) ellei ole kuitua. C Tukkiosan pituus = HPTAB(2) - HPTAB(1) C Kuituosan pituus = HPTAB(3) - HPTAB(2) C Latvahukkaosan pituus = H - HPTAB(3) C HPTAB:ssa olevat mitat ovat IMIT:n mukaan joko maan- C pinnasta tai juurenniskalta. C C Tulosparametrit: C C VTAB(1:5) VTAB(1) Tukkiosan tilavuus, kuutio-dm C VTAB(2) Kuituosan tilavuus, kuutio-dm C VTAB(3) Latvahukkaosan tilavuus, kuutio-dm C VTAB(4) Mahdollisen tyveyksen tilavuus, kuutio-dm. C Olkoon HKAN se arvo joka ohjelma laskee C HPTAB(1):lle ellei sit{ anneta, t.s. kannon C korkeus. Tyveys on HKAN:n ja HPTAB(1):n v{linen C osa edellytt{en ett{ HKAN < HPTAB(1). C VTAB(5) Rungon tilavuus VTAB(1)+VTAB(2)+VTAB(3)+VTAB(4), C t.s. tilavuus korkeudelta min(HKAN,HPTAB(1)) C latvaan. C C Jos sy|tt|parametreiss{ on vikaa on paluutilanteessa jokainen C VTAB(I) = 0.0. HPTAB:ssa on ainakin samat tiedot kuin tulotilan- C teessa, ehk{ muutakin. C C METLA, Matemaattinen osasto - Runkok{yr{t 5.5.1982 / C-G S C FORTRAN IV - VERSIO 2.6.1982 / C-G S C C SUBROUTINE CRKOSA (IPL,IMIT,D,D6,H,MASK,HPTAB,VTAB) C INTEGER IPL, IMIT, MASK REAL D, D6, H, HPTAB(3), VTAB(5) C REAL C(10), DTEOR(3), DKUI(3) C DATA DTEOR /14.5,17.0,16.5/ DATA DKUI /6.3,6.5,6.5/ C C Nollataan tilavuustaulukko ja tarkistetaan ett{ puulaji on sallittu: C DO 1 I = 1,5 VTAB(I) = 0.0 1 CONTINUE C IF (IPL.LT.1 .OR. IPL.GT.3) RETURN C C Estimoidaan juurenniskan korkeus. Mittausten alkukohta on joko C maanpinta tai juurenniska: C IF (IMIT .EQ. 2) GO TO 2 HJN = AKANT1(IPL,D,H) HMIT = 0.0 GO TO 3 2 CONTINUE HJN = HKANIO(IPL,D,H) HMIT = HJN 3 CONTINUE C C Lasketaan runkok{yr{: C CALL CRK(IPL,D,D6,H,HMIT,C) IF (C(9) .EQ. 0.0) RETURN C C Tukkiosan alkukohdaksi asetetaan juurenniska (tai 10 cm maanpinnasta C jos t{m{ on korkeampi) ellei alkukohtaa ole annettu: C HKAN = AMAX1(HJN,0.1) - HMIT IF (MASK .LT. 100) HPTAB(1) = HKAN C C Tukkiosan p{{ttymisl{pimitta on (ellei sit{ ole annettu) m{{r{tt{v{ C PMP-systeemin p|lkytysalgoritmilla. Tarkistetaan kuitenkin aluksi C ett{ puu on tyvess{ niin paksu ett{ yleens{ kannattaa yritt{{ C laskea tukkiosaa: C IF (MOD(MASK,100) .GE. 10) GO TO 4 HPTAB(2) = HPTAB(1) IF (CRKD(HJN,C) .LT. DTEOR(IPL)+1.0) GO TO 4 C C Lasketaan ensin teoreettinen p{{ttymiskorkeus (I-O). Siit{ v{henne- C t{{n EROEST-funktion antama m{{r{ ja korjataan saatua tukkiosapi- C tuutta siten ett{ saadaan s{{nn|nmukaisia p|lkkyj{. Tukkiosa p{{ttyy C kuitenkin viimeist{{n korkeuteen jolla l{pimitta on teoreettinen C p{{ttymisl{pimitta - 1.33 cm: C HTEOR = CRKH(DTEOR(IPL),C) - HJN HPMAA = HJN + HTEOR - EROEST(IPL,HTEOR,NTUK) C IF (CRKD(HPMAA,C) .LT. DTEOR(IPL)-1.33) $ HPMAA = CRKH(DTEOR(IPL)-1.33,C) HT = HPMAA - (HMIT+HPTAB(1)) CALL POLK(HT,HTOD,NTUK) HPTAB(2) = HPTAB(1) + HTOD 4 CONTINUE C C Ellei kuituosan p{{ttymiskorkeutta ole annettu, lasketaan se siten C ett{ puulla on p{{ttymiskorkeudella tietty l{pimitta. Ellei tukki- C osaa saatu on kuituosan kuitenkin oltava v{hint{{n kaksi metri{: C IF (MOD(MASK,10) .NE. 0) GO TO 5 HPTAB(3) = HPTAB(2) IF (CRKD(HJN,C) .LE. DKUI(IPL)) GO TO 5 HKUI = CRKH(DKUI(IPL),C) - (HMIT+HPTAB(2)) IF (HPTAB(2).LE.HPTAB(1) .AND. HKUI.LT.2.0) HKUI = 0.0 HPTAB(3) = HPTAB(2) + HKUI 5 CONTINUE C C Jos HPTAB:ssa ei ole ei-laskeva jono, tilavuuksia ei lasketa: C IF (HPTAB(1).GT.HPTAB(2) .OR. HPTAB(2).GT.HPTAB(3) .OR. $ HPTAB(3).GT.H) RETURN C C Lasketaan osien tilavuudet: C IF (HPTAB(1) .LT. HPTAB(2)) $ VTAB(1) = CRKV(HMIT+HPTAB(1),HMIT+HPTAB(2),C) IF (HPTAB(2) .LT. HPTAB(3)) $ VTAB(2) = CRKV(HMIT+HPTAB(2),HMIT+HPTAB(3),C) VTAB(3) = CRKV(HMIT+HPTAB(3),HMIT+H,C) IF (HKAN .LT. HPTAB(1)) $ VTAB(4) = CRKV(HMIT+HKAN,HMIT+HPTAB(1),C) VTAB(5) = VTAB(1) + VTAB(2) + VTAB(3) + VTAB(4) C RETURN END C ********* C * CRKPK * C ********* C C ALIOHJELMA VALITSEE PERUSRUNKOK[YR[N KERTOIMET PUULAJIN 'IPL' C PERUSTEELLA; IPL = 1, 2 TAI 3. KERTOIMET C(1) - C(8) OVAT C MUUTTUJILLE X, X**2, X**3, X**5, X**8, X**13, X**21 JA X***34. C MATEMAATTINEN OSASTO - RUNKOK[YR[T 27.8.1980/18.5.1981 / C-G S C LEHTIKUUSI LIS[TTY KOODILLA 7 17.9.1982 / C-G S C C SUBROUTINE CRKPK (IPL,C) C DIMENSION C(10) C GO TO (1,2,3,3,5,6,7,1,3), IPL RETURN C C M[NTY: C 1 C(1) = 2.1288 C(2) = -0.63157 C(3) = -1.6082 C(4) = 2.4886 C(5) = -2.4147 C(6) = 2.3619 C(7) = -1.7539 C(8) = 1.0817 RETURN C C KUUSI: C 2 C(1) = 2.3366 C(2) = -3.2684 C(3) = 3.6513 C(4) = -2.2608 C(5) = 0. C(6) = 2.1501 C(7) = -2.7412 C(8) = 1.8876 RETURN C C KOIVU: C 3 C(1) = 0.93838 C(2) = 4.1060 C(3) = -7.8517 C(4) = 7.8993 C(5) = -7.5018 C(6) = 6.3863 C(7) = -4.3918 C(8) = 2.1604 RETURN C C HAAPA: C 5 C(1) = 1.0600 C(2) = 3.7233 C(3) = - 7.4110 C(4) = 7.3395 C(5) = -6.7216 C(6) = 5.9156 C(7) = -4.3415 C(8) = 1.9663 RETURN C C LEPP[: C 6 C(1) = 1.2862 C(2) = 1.4068 C(3) = -2.9250 C(4) = 2.6134 C(5) = -2.5468 C(6) = 2.9720 C(7) = -2.7061 C(8) = 1.5434 RETURN C C LEHTIKUUSI: C 7 C(1) = 1.9046 C(2) = 0.49028 C(3) = -3.4720 C(4) = 4.5721 C(5) = -4.4098 C(6) = 3.7186 C(7) = -2.5256 C(8) = 1.4112 RETURN C END C ******** C * CRKT * C ******** C C OLKOON HX PISTEEN X ET[ISYYS MAANPINNAN TASOSTA (HX:N LAATU: M). C FUNKTIO LASKEE SUHTEEN DX / D80 MISS[ DX ON L[PIMITTA PISTEESS[ X C JA D80 ON L[PIMITTA PISTEESS[, JONKA ET[ISYYS LATVASTA ON 0.8*H. C C(1) - C(8) OVAT RUNKOK[YR[N KERTOIMET JA C(10) PUUN PITUUS H (M). C MATEMAATTINEN OSASTO - RUNKOK[YR[T 27.8.1980 / C-G S C C REAL FUNCTION CRKT (HX,C) C DIMENSION C(10) C X = (C(10)-HX) / C(10) IF (X.GT.1.) GO TO 1 IF (X.LT.0.) GO TO 2 C X2 = X * X X3 = X * X2 X5 = X2 * X3 X8 = X3 * X5 X13 = X5 * X8 X21 = X8 * X13 X34 = X13 * X21 C CRKT = C(1)*X + C(2)*X2 + C(3)*X3 + C(4)*X5 + $ C(5)*X8 + C(6)*X13 + C(7)*X21 + C(8)*X34 RETURN C 1 WRITE (5,601) HX 601 FORMAT (' CRKT: KORKEUS',F8.2,' M < 0') GO TO 3 2 WRITE (5,602) HX, C(10) 602 FORMAT (' CRKT: KORKEUS',F8.2,' M > PUUN PITUUS',F7.2,' M') 3 CRKT = 0. RETURN END C ******** C * CRKV * C ******** C C FUNKTIO LASKEE RUNGON TILAVUUDEN V[LILL[ INTEGROIMALLA C RUNKOK[YR[[. LAADUT: H1 JA H2 (M), CRKV (DM**3). C C(1)-C(8) OVAT RUNKOK[YR[N KERTOIMET, C(9) = D80 ELI PUUN L[PIMITTA C (CM) KORKEUDELLA 0.2*H JA C(10) ON PUUN PITUUS H (M). C MATEMAATTINEN OSASTO - RUNKOK[YR[T 28.8.1980 / C-G S C C REAL FUNCTION CRKV (H1,H2,C) C DIMENSION C(10) DIMENSION A(8,8) DIMENSION HX(2), VOL(2) C H = C(10) C C M[[R[T[[N INTEGROIMISRAJAT HX(1), HX(2) SITEN ETT[ HX(1) < HX(2): C HX(1) = H1 HX(2) = H2 IF (H1 .LT. H2) GO TO 1 HX(1) = H2 HX(2) = H1 C 1 IVIR = 0 IF (HX(1) .LT. 0.) IVIR = 1 IF (HX(2) .GT. H) IVIR = IVIR + 2 IF (IVIR .GT. 0) GO TO 10 C C INTEGROINNISSA TARVITAAN KERTOIMIEN KESKIN[ISET TULOT: C 2 DO 3 I = 1,8 DO 3 J = I,8 A(I,J) = C(I) * C(J) 3 CONTINUE C C LASKETAAN PRIMITIIVIFUNKTION ARVOT VOL(1) JA VOL(2) KORKEUKSIA C HX(1) JA HX(2) VASTAAVISSA RUNKOK[YR[PISTEISS[. LATVASSA (X=0 ELI C HX=H) ON PRIMITIIVIFUNKTION ARVO 0: C DO 4 I = 1,2 X = (H - HX(I)) / H VOL(I) = 0. IF (X .EQ. 0.) GO TO 4 C X2 = X * X X3 = X * X2 X5 = X2 * X3 X8 = X3 * X5 X13 = X5 * X8 C C SEURAAVA LAUSE ON JAETTU KAHTEEN OSAAN SIKSI ETTEI K[[NT[J[ C PYSTY KOKO LAUSETTA K[[NT[M[[N! C VOL(I) = $ ((((( ((((( (((( $ (X13*A(8,8)/69. + A(7,8)/28.) * X8 + A(6,8)/24.) * X5 + $ (2.*A(5,8)+A(7,7))/43.) * X3 + A(4,8)/20.) * X2 + $ A(3,8)/19.) * X + A(2,8)/18.5) * X + $ A(1,8)/18.) * X + A(6,7)/17.5) * X5 + $ A(5,7)/15.) * X3 + $ (2.*A(4,7)+A(6,6))/27.) * X2 + A(3,7)*0.08) * X + $ A(2,7)/12.) * X + A(1,7)/11.5) * X + $ A(5,6)/11.) * X3 + A(4,6)/9.5) * X2 VOL(I) = $ ((((( ((((( (((( VOL(I) + $ (2.*A(3,6)+A(5,5))/17.) * X + A(2,6)*0.125)* X + $ A(1,6)/7.5) * X + A(4,5)/7.) * X2 + $ A(3,5)/6.) * X + $ (2.*A(2,5)+A(4,4))/11.) * X + A(1,5)*0.2) * X + $ A(3,4)/4.5) * X + A(2,4)*0.25) * X + $ (2.*A(1,4)+A(3,3))/7.) * X + A(2,3)/3.) * X + $ (A(1,3)*0.4+A(2,2)*0.2)) * X + A(1,2)*0.5) * X + $ A(1,1)/3.) * X3 4 CONTINUE C C LASKETAAN RUNKO-OSAN TILAVUUS: C V = VOL(1) - VOL(2) CRKV = 0.07853982 * C(9)**2 * H * V RETURN C C VIRHEILMOITUKSET: C 10 IF (IVIR .EQ. 2) GO TO 11 WRITE (5,601) HX(1) 601 FORMAT (' CRKV: KORKEUS',F8.2,' M < 0 - KORVATTU 0:LLA') HX(1) = 0. C 11 IF (IVIR .EQ. 1) GO TO 2 WRITE (5,602) HX(2), H 602 FORMAT (' CRKV: KORKEUS',F8.2,' M KORVATTU PUUN PITUUDELLA', $ F7.2,' M') HX(2) = H GO TO 2 C END C***EROEST C---METLA/MAT/TJP PMPV 1 4.10.1981 FUNCTION EROEST(IPL,HTEOR,NTUK) C---MAARAA KESKIMAARAISEN LATVAKUITUUN SIIRTYVAN TUKKIOSAN C PITUUDEN C---IPL = PUULAJI C---HTEOR = TEOREETTISEN TUKKIOSAN PITUUS C---NTUK = TUKKIEN LUKUMAARAN ESTIMAATTI C--- C---LATVAKUITUUN SIIRTYVAN TUKKIOSAN PITUUS JA TUKKIEN LUKUMAARA C---ESTIMOIDAAN SEURAAVASTI C--- C---SIIRTYMA = AL(I,IPL) + AK(I,IPL)*HTEOR , C---TUKKIEN LKM = I, C--- C---KUN P(I-1,IPL) < HTEOR <= P(I,IPL). C--- C---KERROINTAULUKOT AL, AK JA P REAL AK(5,2), AL(5,2), P(4,2) C--- DATA AK /0.752090,0.689751,0.641760,0.681488,0.652631, +0.848382,0.610232,0.609285,0.617363,0.552112/ DATA AL /3.75490,6.38580,8.79872,12.4718,14.8290, +4.51590,5.74381,8.45665,11.3582,12.4288/ DATA P /8.12758,13.0018,18.1586,24.4122, +7.94745,12.8365,17.3885,22.1969/ C--- C---ALKUASETUKSET ERO=0. NTUK=1 C---FUNKTION ARVO LASKETAAN VAIN PUULAJEILLE IPL=1 JA 2, MUILLE C---EROEST <- 0. IF(IPL.GT.3)GO TO 1000 IF(IPL.EQ.3)GO TO 300 C---MAARAA TUKKIEN LKM DO 100 I=1,4 IF(HTEOR.LE.P(I,IPL))GO TO 200 NTUK=NTUK+1 100 CONTINUE C---ERON ESTIMAATTI SAHALAITAFUNKTIOLTA. C---EROSTA VAHENNETAAN 0.15 M, KOSKA POLKYTYSALGORITMISSA C---TUKKIOSAN C---PITUUS KATKAISTAAN 0.3 M:N MONIKERROIKSI. C---TOIMENPITEEN YHTEISVAIKUTUKSENA TUKKIOSAN PITUUS PYORISTYY C---LAHIMPAAN 0.3 M: MONIKERTAAN. 200 ERO=HTEOR*AK(NTUK,IPL)-AL(NTUK,IPL)-0.15 C---YKSITUKKISILLA RUNGOILLA TUKKIOSAA KASVATETAAN KORKEINATAAN C---1 M. IF(NTUK.LE.1)ERO=AMAX1(ERO,-1.0) GO TO 1000 C---KOIVUILLE KESKIMAARAINEN LATVAKUITUUN SIIRTYVA OSUUS C---ASETETAAN C---1 METRIKSI (L-O 6.7.1981) 300 ERO=1.0 C--- 1000 EROEST=ERO RETURN END C ****** C HKANIO C ****** C FUNKTIO ENNUSTAA JUURENNISKAN KORKEUDEN ILVESSALON MITTOJEN C AVULLA C 8.10.1980/11.10.1982 / C-G S C FUNCTION HKANIO (IPL,D,H) C REAL HKAN(6), DKAN(6) INTEGER INDPUU(9) C DATA HKAN /0.00092874,0.0056579,0.0048152, +0.00011769,0.0017146,0.00805242/ DATA DKAN /0.0045456, 0.0051911,0.0051505, +0.0072661,0.0075402,0.00107911/ DATA INDPUU /1,2,3,3,4,5,6,1,3/ C HKANIO = HKAN(INDPUU(IPL))*H + DKAN(INDPUU(IPL))*D RETURN END C **** C POLK C **** C---METLA/MAT/TJP PMPV 1 4.10.1981 C---SOPIVIA P[TKI[ VALINNUT 5.10.1981 / C-G S SUBROUTINE POLK(PIT,YHTPIT,NPOLK) C---RUNGON OSAN POLKYTTAMINEN POLKYTYSPARAMETRIEN MUKAISESTI. C--- PIT ON P\LKYTT[V[ PITUUS (EROEST:IN J[LKEEN) C--- YHTPIT ON SAATUJEN P\LKKYJEN YHTEISPITUUS C--- NPOLK ON SAATUJEN P\LKKYJEN LUKUM[[R[ C--- C---ST-LAADULLE POLKYTYS ALGORITMI ON SAADETTY C---EROEST-ALIOHJELMALLA C---MAARITETYN YHTENAISEN ST-OSAN PITUUDEN MUKAISEKSI. C---JOS TUKKIOSAN PITUUDEN MAARAAMISTA MUUTETAAN, NIIN TALLOIN ON C---MUUTETTAVA MYOS POLKYTYS ALGORITMIA. C---MUILLE POLKYTETTAVILLE LAADUILLE MAARATAAN POLKKYJEN C---LUKUMAARA C---MINIMIPITUUDEN AVULLA. MAHDOLLINEN YLIJAAMA TASATAAN C---POLKYILLE C---TYVESTA ALKAEN KASVATTAMALLA TARPEELLINEN MAARA POLKKYJA C---MAKSIMIPITUISIKSI. C--- C---ASETTAA C---OSITUS/NPOLK <- SAATUJEN POLKKYJEN LUKUMAARA C--- /POLPIT (V) <- SAATUJEN POLKKYJEN PITUUDET C--- C---KUTSUTUT ALIOHJELMAT C--- C--- DIMENSION POLPIT(10) C PMIN = 4.3 PKES = 4.9 PMAX = 6.1 C C--- C---PLK1: ALKUASETUKSET, ST LAADUILLE 0-3 POLKYN PITUUDEN TULEE OLLA C---0.3 M:N MONIKERTOJA; MUILLA LAADUILLA 0.1 M:N MONIKERTOJA. C---JAKOVALI ASETETAAN TAMAN MUKAISESTI. NPOLK=0 YHTPIT = 0 VALI=0.3 C---PLK2: TARKISTA TULEEKO MINIMIPOLKKY, MIKALI ST-LAADUSTA EI TULE C---PARAMETRINA ANNETTUA MINIMIPOLKKYA, NIIN SIITA YRITETAAN C---LYHYTTA C---POLKKYA YLEISTEN LAPIMITTAVAATIMUSTEN MUKAISESTI. IF(PIT.GE.PMIN)GO TO 300 GO TO 9000 C---PLK3: YHDEN POLKYN NORMAALIT RUNGON OSAT 300 IF(PIT.GT.PMAX)GO TO 400 NPOLK=1 POLPIT(1)=PIT-AMOD(PIT-0.099,VALI)+0.001 YHTPIT=POLPIT(1) GO TO 9000 C---PLK4: USEAN POLKYN RUNGON OSAT. C---MAARAA POLKKYJEN LKM SEKA POLKYN PITUUS. 400 NPOLK=PIT/PKES+0.5 IF(PIT/NPOLK.LT.PMIN)NPOLK=NPOLK-1 PPIT=AMIN1(PIT/NPOLK,PMAX) PPIT=PPIT-AMOD(PPIT-0.099,VALI)+0.001 C---PLK5: TASAA YLIJAAMA TYVESTA ALKAEN KASVATTAMALLA POLKYT C---MAKSIMIPITUISIKSI, JOS TARPEEN ERO=PIT-NPOLK*PPIT ERO=ERO-AMOD(ERO+0.001,VALI)+0.001 DO 510 I=1,NPOLK DEL=0. IF(ERO.GT.0)DEL=AMIN1(PMAX-PPIT,ERO) POLPIT(I)=PPIT+DEL ERO=ERO-DEL YHTPIT=YHTPIT+POLPIT(I) 510 CONTINUE C---LOPETA 9000 RETURN END C--- C ******** C * CV2K * C ******** C C FUNKTIO LASKEE PUUN KUORELLISEN TILAVUUDEN (KUUTIO-DM) MAASTA C MITATTUJEN TUNNUSTEN C D RINNANKORKEUSL[PIMITTA (CM), 0.1 <= D <= 85. JA C H PUUN PITUUS (M), 1.4 <= H <= 42. C PERUSTEELLA. C PUULAJI ON M[NTY (IPL=1), KUUSI (IPL=2) TAI KOIVU (IPL=3 TAI 4). C LEHTIKUUSI (IPL=7) LIS[TTY 4.3.1983 / C-G S C LEHTIKUUSEN KOODI MUUTETTU 7 -> 9 13.5.1984 / C-G S C MATEMAATTINEN OSASTO - TILAVUUSFUNKTIOT 13.10.1980/4.3.1981 / C-G S C Pienille puille lis{tty Antti Ihalaisen laskemat erikoisfunktiot C 5.2.1986 / C-G S C C FUNCTION CV2K (IPL,D,H) C REAL HRAJA(4), C1(3), C2(3), C3(3), C4(3) LOGICAL V1, V2, V3 C DATA HRAJA /4.5,3.5,6.5,6.5/ DATA C1 /-4.21324,-3.44416,-4.92615/ DATA C2 /2.33888,1.78061,2.56275/ DATA C3 /0.675080,0.708673,0.787837/ DATA C4 /0.0100470,0.0740863,-0.0106150/ C C TUNNUSTEN TARKISTUS: C V1 = (IPL.LT.1 .OR. IPL.GT.9) .AND. IPL.NE.9 V2 = D.LT.0.1 .OR. D.GT.85. V3 = H.LT.1.4 .OR. H.GT.42. IF (V1 .OR. V2 .OR. V3) GO TO 20 C C PIENET PUUT (IHALAINEN): C IF (IPL .LE. 4) THEN IF (H .LE. HRAJA(IPL)) THEN HLN = ALOG(H) DXLN = ALOG(2.0+1.25*D) IP = MIN(IPL,3) VLN = C1(IP) + C2(IP)*DXLN + C3(IP)*HLN + C4(IP)*D GO TO 19 END IF END IF C C APUMUUTTUJAT: C DLN = ALOG(D) HLN = ALOG(H) H13LN = ALOG(H-1.3) C C TILAVUUDET PUULAJEITTAIN: C GO TO (11,12,13,13,13,13,13,12,17), IPL C C M[NTY: C 11 VLN = - 3.32176 + 2.01395*DLN + 2.07025*HLN - 1.07209*H13LN $ - 0.0032473*D GO TO 19 C C KUUSI: C 12 VLN = - 3.77543 + 1.91505*DLN + 2.82541*HLN - 1.53547*H13LN $ - 0.0085726*D GO TO 19 C C KOIVU: C 13 VLN = - 4.49213 + 2.10253*DLN + 3.98519*HLN - 2.65900*H13LN $ - 0.014097*D GO TO 19 C C LEHTIKUUSI: C 17 VLN = - 4.30421 + 2.06818*DLN + 3.53310*HLN - 2.23567*H13LN $ - 0.0111362*D C 19 CV2K = EXP(VLN) RETURN C C VIRHEILMOITUKSET: C 20 WRITE (6,600) IF (V1) WRITE (6,601) IPL IF (V2) WRITE (6,602) D IF (V3) WRITE (6,603) H CV2K = 0.0 RETURN C 600 FORMAT (' CV2K: VIRHEELLISET ARGUMENTIT!') 601 FORMAT (9X,'PUULAJITUNNUS =',I3,' 1, 2, 3, 4 JA 9 SALLITAAN') 602 FORMAT (9X,'L[PIMITAN D (',F8.2,') ON OLTAVA V[LILL[ (0.1,85.)') 603 FORMAT (9X,'PITUUDEN H (',F8.2,') ON OLTAVA V[LILL[ (1.4,42.)') C END C ******** C * CV3K * C ******** C C FUNKTIO LASKEE PUUN KUORELLISEN TILAVUUDEN (KUUTIO-DM) MAASTA C MITATTUJEN TUNNUSTEN C D RINNANKORKEUSL[PIMITTA (CM), 0.1 <= D <= 85., C D6 KUUDEN METRIN L[PIMITTA (CM), 0.1 <= D6 <= D JA C H PUUN PITUUS (M), 6.1 <= H <= 40. C PERUSTEELLA. C PUULAJI ON M[NTY (IPL=1), KUUSI (IPL=2) TAI KOIVU (IPL=3 TAI 4). C KOIVUN KOODI 4 LIS[TTY 13.5.1984 / C-G S C LEHTIKUUSI, IPL=9, LIS[TTY 24.2.1987 / C-G S C MATEMAATTINEN OSASTO - TILAVUUSFUNKTIOT 13.10.1980/4.3.1981 / C-G S C C FUNCTION CV3K (IPL,D,D6,H) C LOGICAL V1, V2, V3, V4 C C TUNNUSTEN TARKISTUS: C V1 = (IPL.LT.1 .OR. IPL.GT.9) .AND. IPL.NE.9 V2 = D.LT.0.1 .OR. D.GT.85. .AND. IPL.NE.9 V3 = D6.LT.0.1 .OR. D6.GT.D V4 = H.LT.6.1 .OR. H.GT.40. .AND. IPL.NE.9 IF (V1 .OR. V2 .OR. V3 .OR. V4) GO TO 20 C C APUMUUTTUJAT: C D2 = D**2 D2H = D2*H D3H = D*D2H D2H2 = D2H*H D62 = D6**2 DSU = D2 + D*D6 + D62 D62H6 = D62 * (H-6.) C C TILAVUUDET PUULAJEITTAIN: C GO TO (11,12,13,13,13,13,13,12,19), IPL C C M[NTY: C 11 CV3K = 0.268621*D2 - 0.0145543*D2H - 0.0000478628*D3H + $ 0.000334101*D2H2 + 0.0973148*DSU + 0.0440716*D62H6 RETURN C C KUUSI: C 12 CV3K = 0.208043*D2 - 0.0149567*D2H - 0.000114406*D3H + $ 0.000436781*D2H2 + 0.133947*DSU + 0.0374599*D62H6 RETURN C C KOIVU: C 13 CV3K = 0.226547*D2 - 0.0104691*D2H - 0.000122258*D3H + $ 0.000438033*D2H2 + 0.0991620*DSU + 0.0334836*D62H6 RETURN C C LEHTIKUUSI: C 19 CV3K = 0.281982*D2 - 0.0152185*D2H - 0.0000980420*D3H + $ 0.000370942*D2H2 + 0.0914879*DSU + 0.0441584*D62H6 RETURN C C VIRHEILMOITUKSET: C 20 WRITE (6,600) IF (V1) WRITE (6,601) IPL IF (V2) WRITE (6,602) D IF (V3) WRITE (6,603) D6, D IF (V4) WRITE (6,604) H CV3K = 0.0 RETURN C 600 FORMAT (' CV3K: VIRHEELLISET ARGUMENTIT!') 601 FORMAT (9X,'PUULAJITUNNUS =',I3,' 1, 2,3 JA 4 SALLITAAN') 602 FORMAT (9X,'L[PIMITAN D (',F8.2,') ON OLTAVA V[LILL[ (0.1,85.)') 603 FORMAT (9X,'L[PIMITAN D6 (',F8.2,') ON OLTAVA V[LILL[ (0.1,D) ', $ 'D =',F8.2) 604 FORMAT (9X,'PITUUDEN H (',F8.2,') ON OLTAVA V[LILL[ (6.1,40.)') C END SUBROUTINE LLASKU(LEIMA,LKML,MALLI,LKMM,D,D6,PIT,TIL,I,PPPA, +PTIL,PDGM,PHGM,JPPA,JTIL,JDGM,JHGM,EPPA,ETIL,EDGM,EHGM,PALA) PARAMETER (MAXP=500) INTEGER LEIMA(MAXP),LKML,MALLI(MAXP),LKMM,I REAL D(MAXP),D6(MAXP),PIT(MAXP),TIL(MAXP),PPPA,PTIL,PDGM,PHGM, +JPPA,JTIL,JDGM,JHGM,PALA,D2,D3,SUMD2,SUMD3,SUMTIL,D2H,SUMD2H, +EPPA,ETIL,EDGM,EHGM C Lasketaan poistuman tunnukset, poistettavien puiden numerot ovat C taulukossa LEIMA, ne ovat indeksej„ puutunnustaulukoille D,D6,PIT ja TIL C D(cm),D6(cm),PIT(m),TIL(dm3),PALA(ha) IF(LKML.EQ.0) GOTO 801 DO 1 IA=1,LKML D2=(D(LEIMA(IA))/100.)**2 D3=(D(LEIMA(IA))/100.)**3 D2H=D2*PIT(LEIMA(IA)) SUMD2H=SUMD2H+D2H SUMTIL=SUMTIL+TIL(LEIMA(IA)) SUMD2=SUMD2+D2 SUMD3=SUMD3+D3 1 CONTINUE C Poistuman pohjapinta-ala PPPA hehtaarikohtaisena A1=(3.14159/4.)*SUMD2 A2=1./PALA PPPA=A1*A2 C Poistuman keskil„pimitta PDGM=SUMD3/SUMD2 C Poistuman keskipituus PHGM=SUMD2H/SUMD2 C Poistuman tilavuus kuutiometri„ hehtaarilla PTIL=(SUMTIL/1000.)*A2 C ************************************************************ C J„„v„n puuston tunnusten laskemiseksi tulisi hakea ne alkiot C taulukoista D D6 PIT ja TIL, jotka eiv„t ole leimattuja C (niihin ei viitata LEIMA-taulukossa) 801 CONTINUE SUMD2H=0. SUMD2=0. SUMD3=0. SUMTIL=0. DO 5 IA=1,I DO 6 IB=1,LKML IF(LEIMA(IB).EQ.IA) GOTO 5 6 CONTINUE D2=(D(IA)/100.)**2 D3=(D(IA)/100.)**3 D2H=D2*PIT(IA) SUMD2H=SUMD2H+D2H SUMTIL=SUMTIL+TIL(IA) SUMD2=SUMD2+D2 SUMD3=SUMD3+D3 5 CONTINUE C J„„v„n puuston pohjapinta-ala JPPA hehtaarikohtaisena A1=0.7853*SUMD2 A2=1./PALA JPPA=A1*A2 C J„„v„n puuston keskilApimitta JDGM=SUMD3/SUMD2 C J„„v„n puuston keskipituus JHGM=SUMD2H/SUMD2 C j„„v„n puuston tilavuus kuutiometri„ hehtaarilla JTIL=(SUMTIL*A2)/1000. C ***************************************************************** C Metsikk”tunnukset ennen harvennusta saadaan k„ym„ll„ l„pi C kaikki I puuta SUMD2H=0. SUMD2=0. SUMD3=0. SUMTIL=0. DO 8 IA=1,I D2=(D(IA)/100.)**2 D3=(D(IA)/100.)**3 D2H=D2*PIT(IA) SUMD2H=SUMD2H+D2H SUMTIL=SUMTIL+TIL(IA) SUMD2=SUMD2+D2 SUMD3=SUMD3+D3 8 CONTINUE C Ennen harvennusta pohjapinta-ala EPPA hehtaarikohtaisena A1=0.7853*SUMD2 A2=1./PALA EPPA=A1*A2 C Ennen harvennusta puuston keskil„pimitta EDGM=SUMD3/SUMD2 C Ennen harvennusta puuston keskipituus EHGM=SUMD2H/SUMD2 C Ennen harvennusta puuston tilavuus kuutiometri„ hehtaarilla ETIL=(SUMTIL*A2)/1000. RETURN END SUBROUTINE PITK(IALKU,JNRO,D,H,B0,B1) C Aliohjelmalla lasketaan arvot parametreille B0 ja B1 pituusk„yr„lle C D:n funktiona. Parametrit ratkaistaan laskemalla useita ratkaisuja C ja valitsemalla niist„ paras. C H=B0*[D/(D+B1)] C C Aliohjelmalle v„litett„v„t parametrit: C ALKU = indeksi, jolla puutietotaulukoiden luku alkaa (D(IALKU),H(IALKU)) C JNRO = indeksi, viimeinen alkio taulukoissa C D = rinnankorkeusl„pimitat (cm) C H = pituudet (m) C Aliohjelma palauttaa parametrit: C B0 = pituusmallin parametri (ratkaistaan) C B1 = pituusmallin parametri (ratkaistaan) C INTEGER I,J,K,L,IX,JX,MAXP,IALKU,JNRO PARAMETER (MAXP=500) REAL ALPHA(20), D(MAXP), H(MAXP), HPRED(MAXP), SUMRES(20,30), + BETA(30) , B0, B1, MINRES, RESID(MAXP) DO 4 I=1,20 ALPHA(I) = 5.00 + REAL(I) * 3.0 4 CONTINUE DO 5 I=1,30 BETA(I) = 0.01 + REAL(I) * 2.0 5 CONTINUE C Tehd„„n 600 kappaletta pituusk„yri„, joiden joukosta valitaan C PNS-estimaattori. (20 x 30 = 600) WRITE (6,*) 'Creating 600 height curves and fitting them.' DO 6 I=1,20 DO 7 J=1,30 SUMRES(I,J)=0. DO 8 L=IALKU,JNRO C I siirt„„ alpha (B0) parametria 3.0 metrin askeleissa. C J siirt„„ beta (B1) parametria 2.0 metrin askeleissa. C L k„y l„pi kaikki koepuut kullakin kuudellasadasta parametrikombinaatiosta. C HPRED -taulukkoon lasketaan puulle L parametrien ALPHA(I) ja BETA(J) C arvoilla k„yr„n ja todellisen pituuden erotus. C RESID -taulukkoon lasketaan puulle L vastaavan erotuksen neli”. C SUMRES -taulukkoon lasketaan yhteen residuaalien summaa parametreille C ALPHA(I) ja BETA(J). HPRED(L) = ALPHA(I) * (D(L) / (D(L) + BETA(J) ) ) RESID(L) = ( H(L) - HPRED(L) ) ** 2 SUMRES(I,J) = SUMRES(I,J) + RESID(L) 8 CONTINUE 7 CONTINUE 6 CONTINUE C Haetaan PNS-estimaattori taulukon SUMRES avulla. ISO=0 WRITE (6,*) ' Residuals computed - wait ! ' MINRES = 9999999999. C Haetaan SUMRES-taulukon pienin alkio, johon viittaavat indeksit C II ja JJ. DO 10 I=1,20 DO 11 J=1,30 IF( SUMRES(I,J) .LT. MINRES) THEN MINRES = SUMRES (I,J) II = I JJ = J END IF 11 CONTINUE 10 CONTINUE B0 = ALPHA(II) B1 = BETA(JJ) DO 16 I=1,20 DO 17 J=1,30 SUMRES(I,J)=0. 17 CONTINUE 16 CONTINUE RETURN END SUBROUTINE REGD6(D,D6,I,A,BE) PARAMETER (MAXP=500) REAL D(MAXP),D6(MAXP),A,BE,KAD,KAD6,SUMD,SUMD6,APUNIM,APUOS, +SUMNIM,SUMOS INTEGER I,IM C D-taulukko v„litt„„ rinnankorkeusl„pimitat C D6-taulukko v„litt„„ yl„l„pimitat C I laskuri kertoo puiden m„„r„n koealalla C A palauttaa regressioyht„l”n vakion C BE palauttaa regressiokertoimen C Lasketaan lineaarinen regressio d1.3 ja d6 v„lille æ = A + BE*X C Haetaan koepuut, joista d6 mitattu ja lasketaan keskiarvot DO 15 M=1,I IF (D6(M).GT.0.1) THEN IM=IM+1 SUMD6=SUMD6+D6(M) SUMD=SUMD+D(M) END IF 15 CONTINUE KAD6=SUMD6/REAL(IM) KAD=SUMD/REAL(IM) DO 16 M=1,I IF (D6(M).GT.0.1) THEN APUNIM=(D(M)-KAD)*(D6(M)-KAD6) SUMNIM=SUMNIM+APUNIM APUOS =(D(M)-KAD)**2. SUMOS=SUMOS+APUOS END IF 16 CONTINUE BE=(SUMNIM/SUMOS) A=KAD6-(BE*KAD) RETURN END SUBROUTINE RESULT(LEIMA,LKML,MALLI,LKMM,I,VP,VJ,TULOS,APU2) PARAMETER (MAXP=500) INTEGER LEIMA(MAXP),LKML,MALLI(MAXP),LKMM,I,IX,JX REAL VP,VJ,TULOS CHARACTER APU2*23,NIMI*23 DO 10 L=12,20 IF (APU2(L:L).EQ.' ') GOTO 9 10 CONTINUE 9 CONTINUE DO 8 M=12,23 IF(APU2(M:M).EQ.'.' .AND. M.LT.L) GOTO 15 8 CONTINUE NIMI(1:23)='C:\LEIMAUS\'//APU2(12:L-1)//'.VIR' GOTO 12 15 CONTINUE NIMI(1:23)='C:\LEIMAUS\'//APU2(12:M-1)//'.VIR' GOTO 12 12 CONTINUE OPEN (17,FILE=NIMI,STATUS='NEW') C I=Puiden kokonaism„„r„, LKMM=mallileimattujen puiden lukum„„r„ C LKML=leimaajan poistamien lukum„„r„ C Lasketaan kilpailijan v„„rin poistamien puiden lukum„„r„ VP=0 DO 1 N=1,LKML IX=0 DO 2 M=1,LKMM IF(LEIMA(N).EQ.MALLI(M)) IX=IX+1 2 CONTINUE IF(IX.EQ.0) THEN VP=VP+1 WRITE (17,201) LEIMA(N) END IF 201 FORMAT (' ',I3,2X,'P') 1 CONTINUE C Lasketaan kilpailijan v„„rin j„tettyjen puiden lukum„„r„ VJ=0 DO 3 N=1,LKMM JX=0 DO 4 M=1,LKML IF(MALLI(N).EQ.LEIMA(M)) JX=JX+1 4 CONTINUE IF(JX.EQ.0) THEN VJ=VJ+1 WRITE (17,202) MALLI(N) 202 FORMAT (' ',I3,2X,'J') END IF 3 CONTINUE C Tulos on puiden kokonaism„„r„ I v„hennettyn„ j„tt”- ja poistovirheiden C lukum„„r„ll„ VJ+VP TULOS=I-(VP+VJ) CLOSE (17) RETURN END SUBROUTINE TEXT +(VPLAJI,VOSUUS,SUMTUK,SUMKUI,SUMHUK,SUMTOT,SUMTIL,RLUKU,RTUK, +RKUI,RPOI,TUKKES,PPPA,PDGM,PTIL,PHGM,JPPA,JDGM,JTIL,JHGM,EPPA, +EDGM,ETIL,EHGM,PALA,D,D6,PIT,TIL,VTUK,VKUI,VHUK,VTOT,I,TULOS, +VP,VJ,LKML,HDOM,JHDOM) include 'fgraph.fd' PARAMETER (MAXP=500) REAL VPLAJI(4),VOSUUS(4),SUMTUK,SUMKUI,SUMHUK,SUMTOT,SUMTIL, +RLUKU,RTUK,RKUI,RPOI,TUKKES,PPPA,PDGM,PTIL,PHGM,JPPA,JDGM,JTIL, +JHGM,EPPA,EDGM,ETIL,EHGM,PALA,D(MAXP),D6(MAXP),PIT(MAXP), +TIL(MAXP),VTUK(MAXP),VKUI(MAXP),VHUK(MAXP),VTOT(MAXP),TULOS, +VP,VJ,A1,RJAA,RLEI,HDOM,JHDOM INTEGER I,LKML CHARACTER RIVI*80,APU*80,NIMI*80,APU1*80 record /rccoord/ curpos 801 CONTINUE CALL clearscreen ($GCLEARSCREEN) CALL settextposition (4,1,curpos) CALL outtext (' Once again - please enter...') CALL settextposition (5,1,curpos) CALL outtext ( ' Your full name (first name and family name:') READ (5,'(A80)') NIMI DO 5 IY=1,80 IF (NIMI(1:5).EQ.' ') GOTO 801 IF (NIMI(IY:IY).EQ.' ' .AND. NIMI(IY:IY+5).EQ.' ') GOTO 1 5 CONTINUE 1 CONTINUE CALL clearscreen ($GCLEARSCREEN) CALL settextposition (4,1,curpos) CALL outtext ('Your results will be written to file defined here') CALL settextposition (5,1,curpos) CALL outtext (' Please enter a shortened version (max. 8 char.)') CALL settextposition (6,1,curpos) CALL outtext (' of your family name (surname):') READ (5,'(A80)') APU1 DO 2 IX=1,80 IF (APU1(IX:IX).EQ.' ') GOTO 3 2 CONTINUE 3 CONTINUE IF (IX.EQ.1) GOTO 1 IF (IX.GT.9) THEN WRITE (6,*) ' 8 characters is the maximum, try again.' GOTO 1 END IF APU(1:IX+14)='C:\LEIMAUS\'//APU1(1:IX-1)//'.TUL' C Leimattuja runkoa LKML kappaletta koealalla C Hehtaarilla RLEI kappaletta ja j„„v„ runkoluku RJAA kappaletta A1=1./PALA RJAA=(REAL(I)-REAL(LKML))*A1 RLEI=REAL(LKML)*A1 OPEN (18,FILE=APU,STATUS='NEW') RIVI(1:30+IY)=' Results of the Competition '//NIMI(1:IY) WRITE(18,*) ' ' WRITE(18,'(A80)') RIVI WRITE(18,*) ' ' WRITE(18,*) ' Variable Value' WRITE(18,*) ' ' WRITE(18,'(A33,2X,F10.2)') ' Stems [r/ha] ............ ', +RLUKU WRITE(18,'(A33,2X,F10.2)') ' Sawlog stems [r/ha] .......... ', +RTUK WRITE(18,'(A33,2X,F10.2)') ' Pulpwood stems [r/ha] ........ ', +RKUI WRITE(18,'(A33,2X,F10.2)') ' IGNORE ...... ', +RPOI WRITE(18,*) ' ' WRITE(18,'(A33,2X,F10.2)') ' Mean volume of saw-log [dm3] ', +TUKKES WRITE(18,'(A33,2X,F10.2)') ' Tot. volume. (CV3K) [m3/ha].. ', +SUMTIL/1000. WRITE(18,'(A33,2X,F10.2)') ' Tot. volume (CRK) [m3/ha] .. ', +SUMTOT/1000. WRITE(18,'(A33,2X,F10.2)') ' Saw-wood [m3/ha] ............ ', +SUMTUK/1000. WRITE(18,'(A33,2X,F10.2)') ' Pulpwood [m3/ha] ............ ', +SUMKUI/1000. WRITE(18,'(A33,2X,F10.2)') ' Non-commercial [m3/ha] ...... ', +SUMHUK/1000. WRITE(18,'(A33,2X,F10.2)') ' Pine [m3/ha] ............... ', +VPLAJI(1)/1000. WRITE(18,'(A33,2X,F10.2)') ' Spruce [m3/ha] ............... ', +VPLAJI(2)/1000. WRITE(18,'(A33,2X,F10.2)') ' Birch [m3/ha] ............... ', +VPLAJI(3)/1000. WRITE(18,'(A33,2X,F10.2)') ' ............... ', +VPLAJI(4)/1000. WRITE(18,'(A33,2X,F10.2)') ' Pine [%] ................... ', +VOSUUS(1) WRITE(18,'(A33,2X,F10.2)') ' Spruce [%] ................... ', +VOSUUS(2) WRITE(18,'(A33,2X,F10.2)') ' Birch [%] ................... ', +VOSUUS(3) WRITE(18,'(A33,2X,F10.2)') ' ................... ', +VOSUUS(4) WRITE(18,*) ' ' WRITE(18,'(A42)') ' Var Before Cut After' WRITE(18,'(A16,F7.2,6X,F7.2,6X,F7.2)')' G [m2/ha] ',EPPA, +PPPA,JPPA WRITE(18,'(A16,F7.2,6X,F7.2,6X,F7.2)')' Dgm [cm] ', +EDGM*100.,PDGM*100.,JDGM*100. WRITE(18,'(A16,F7.2,6X,F7.2,6X,F7.2)')' V [m3/ha] ',ETIL, +PTIL,JTIL WRITE(18,'(A16,F7.2,6X,F7.2,6X,F7.2)')' Hgm [cm] ',EHGM, +PHGM,JHGM WRITE(18,'(A16,F7.2,6X,F7.2,6X,F7.2)')' Stems [r/ha] ',RLUKU, +RLEI,RJAA WRITE(18,'(A16,F7.2,6X,13X,F7.2)') ' Hdom [m] ',HDOM, +JHDOM WRITE (18,*) WRITE (18,'(A40,F4.0)') ' Your RESULT is ',TULOS WRITE (18,'(A40,F4.0)') ' Should not have left ',VJ WRITE (18,'(A40,F4.0)') ' Should not have been removed ',VP WRITE (18,'(A40,F4.0)') ' Max result could have been ',REAL(I) CLOSE (18) CALL clearscreen ($GCLEARSCREEN) CALL settextposition (4,1,curpos) CALL outtext('Your results have been written to ') CALL settextposition (5,1,curpos) CALL outtext(APU(1:25)) CALL settextposition (6,1,curpos) CALL outtext(' file') CALL settextposition (7,1,curpos) CALL outtext('in the C:\LEIMAUS -subdirectory.') CALL settextposition (9,1,curpos) CALL outtext(' Complete the program and return to DOS.') CALL settextposition (10,1,curpos) CALL outtext(' Results are printed to the HP-laser printer using') CALL settextposition (11,1,curpos) CALL outtext(' the command:') status = settextcolor (6) CALL outtext(' PSTXT ') CALL outtext(APU(1:26)) CALL outtext(' MIKROLUOKKA') CALL settextposition (14,1,curpos) CALL outtext(' Press return to continue...') READ (*,*) RETURN END SUBROUTINE TILOS(IPL,VTOT,I,VPLAJI,VOSUUS,PALA) PARAMETER (MAXP=500) REAL VTOT(MAXP),VPLAJI(4),VOSUUS(4),A1,SUMV INTEGER IPL(MAXP),I,IIPL,J C VTOT(MAXP) runkok„yrill„ lasketut runkotilavuudet C I puiden lukum„„r„ koealalla C IPL puulajitaulukko (Huom. CRKOSA laskenut muulp koivuna!) C VOSUUS laskentaryhmien tilavuusosuudet prosentteina DO 1 J=1,I IIPL=IPL(J) IF (IPL(J).GE.3 .AND. IPL(J).LE.4) IIPL=3 IF (IPL(J).GE.5 .AND. IPL(J).LE.7) IIPL=4 IF (IPL(J).EQ.8) IIPL=2 IF(IIPL.EQ.1) VPLAJI(1)=VPLAJI(1)+VTOT(J) IF(IIPL.EQ.2) VPLAJI(2)=VPLAJI(2)+VTOT(J) IF(IIPL.EQ.3) VPLAJI(3)=VPLAJI(3)+VTOT(J) IF(IIPL.EQ.4) VPLAJI(4)=VPLAJI(4)+VTOT(J) SUMV=SUMV+VTOT(J) 1 CONTINUE A1=1./PALA C Muunnetaan tulokset hehtaarikohtaisiksi ja lasketaan puulajiosuudet C SUMV kokonaistilavuus koealalla SUMV=SUMV*A1 DO 2 J=1,4 VPLAJI(J)=VPLAJI(J)*A1 VOSUUS(J)=(VPLAJI(J)/SUMV)*100. 2 CONTINUE RETURN END C Aliohjelma KILP-ohjelmaan, joka tulostaa C:\HARSI\ hakemistoon C HARSIa varten PUUSTO.DAT tiedoston PUUSTO.HAR-nimisen„ ja oikeassa C formaatissa sek„ *.LEI leimattujen puiden numerot harsia varten SUBROUTINE TTEKO(NRO,X,Y,IPL,D,PIT,I,APU2,LEIMA,LKML) PARAMETER (MAXP=500) INTEGER NRO(MAXP),X(MAXP),Y(MAXP),IPL(MAXP),I,IIPL(MAXP),DI(MAXP), +LEIMA(MAXP),LKML,IXMAX,IYMAX,IYMIN,IXMIN REAL D(MAXP),PIT(MAXP),PIIT(MAXP) CHARACTER APU2*23,ULEI*26,UHAR*23 LOGICAL LINQ DO 3 K=12,23 IF(APU2(K:K).EQ.'.') GOTO 4 IF(APU2(K:K).EQ.' ') GOTO 4 3 CONTINUE 4 CONTINUE ULEI(1:26)='C:\LEIMAUS\'//APU2(12:K-1)//'.LEI' OPEN (20,FILE=ULEI,STATUS='NEW') C Pituuksien py”rist„minen l„himp„„n 0.5 metriin if-lauseilla C L„pimittojen py”rist„minen parittomiin 2 cm:n luokkiin C PIT -> PIIT aliohjelma muuttaa pituudet DO 1 J=1,I IF(PIT(J).LT.1.0) PIIT(J)=0.5 IF(PIT(J).GE.1.0 .AND. PIT(J).LT.2.0) PIIT(J)=1.5 IF(PIT(J).GE.2.0 .AND. PIT(J).LT.3.0) PIIT(J)=2.5 IF(PIT(J).GE.3.0 .AND. PIT(J).LT.4.0) PIIT(J)=3.5 IF(PIT(J).GE.4.0 .AND. PIT(J).LT.5.0) PIIT(J)=4.5 IF(PIT(J).GE.5.0 .AND. PIT(J).LT.6.0) PIIT(J)=5.5 IF(PIT(J).GE.6.0 .AND. PIT(J).LT.7.0) PIIT(J)=6.5 IF(PIT(J).GE.7.0 .AND. PIT(J).LT.8.0) PIIT(J)=7.5 IF(PIT(J).GE.8.0 .AND. PIT(J).LT.9.0) PIIT(J)=8.5 IF(PIT(J).GE.9.0 .AND. PIT(J).LT.10.0) PIIT(J)=9.5 IF(PIT(J).GE.10.0 .AND. PIT(J).LT.11.0) PIIT(J)=10.5 IF(PIT(J).GE.11.0 .AND. PIT(J).LT.12.0) PIIT(J)=11.5 IF(PIT(J).GE.12.0 .AND. PIT(J).LT.13.0) PIIT(J)=12.5 IF(PIT(J).GE.13.0 .AND. PIT(J).LT.14.0) PIIT(J)=13.5 IF(PIT(J).GE.14.0 .AND. PIT(J).LT.15.0) PIIT(J)=14.5 IF(PIT(J).GE.15.0 .AND. PIT(J).LT.16.0) PIIT(J)=15.5 IF(PIT(J).GE.16.0 .AND. PIT(J).LT.17.0) PIIT(J)=16.5 IF(PIT(J).GE.17.0 .AND. PIT(J).LT.18.0) PIIT(J)=17.5 IF(PIT(J).GE.18.0 .AND. PIT(J).LT.19.0) PIIT(J)=18.5 IF(PIT(J).GE.19.0 .AND. PIT(J).LT.20.0) PIIT(J)=19.5 IF(PIT(J).GE.20.0 .AND. PIT(J).LT.21.0) PIIT(J)=20.5 IF(PIT(J).GE.21.0 .AND. PIT(J).LT.22.0) PIIT(J)=21.5 IF(PIT(J).GE.22.0 .AND. PIT(J).LT.23.0) PIIT(J)=22.5 IF(PIT(J).GE.23.0 .AND. PIT(J).LT.24.0) PIIT(J)=23.5 IF(PIT(J).GE.24.0 .AND. PIT(J).LT.25.0) PIIT(J)=24.5 IF(PIT(J).GE.25.0 .AND. PIT(J).LT.26.0) PIIT(J)=25.5 IF(PIT(J).GE.26.0 .AND. PIT(J).LT.27.0) PIIT(J)=26.5 IF(PIT(J).GE.27.0 .AND. PIT(J).LT.28.0) PIIT(J)=27.5 IF(PIT(J).GE.28.0 .AND. PIT(J).LT.29.0) PIIT(J)=28.5 IF(PIT(J).GE.29.0 .AND. PIT(J).LT.30.0) PIIT(J)=29.5 IF(PIT(J).GE.30.0 .AND. PIT(J).LT.31.0) PIIT(J)=30.5 IF(PIT(J).GE.31.0 .AND. PIT(J).LT.32.0) PIIT(J)=31.5 IF(PIT(J).GE.32.0 .AND. PIT(J).LT.33.0) PIIT(J)=32.5 IF(PIT(J).GE.33.0 .AND. PIT(J).LT.34.0) PIIT(J)=33.5 IF(PIT(J).GE.34.0 .AND. PIT(J).LT.35.0) PIIT(J)=34.5 IF(PIT(J).GE.35.0 .AND. PIT(J).LT.36.0) PIIT(J)=35.5 IF(PIT(J).GE.36.0 .AND. PIT(J).LT.37.0) PIIT(J)=36.5 IF(PIT(J).GE.37.0 .AND. PIT(J).LT.38.0) PIIT(J)=37.5 IF(PIT(J).GT.38.0) PIIT(J)=38.5 IF(D(J).LT.2.0) DI(J)=1 IF(D(J).GE.2.0 .AND. D(J).LT.4.0) DI(J)=3 IF(D(J).GE.4.0 .AND. D(J).LT.6.0) DI(J)=5 IF(D(J).GE.6.0 .AND. D(J).LT.8.0) DI(J)=7 IF(D(J).GE.8.0 .AND. D(J).LT.10.0) DI(J)=9 IF(D(J).GE.10.0 .AND. D(J).LT.12.0) DI(J)=11 IF(D(J).GE.12.0 .AND. D(J).LT.14.0) DI(J)=13 IF(D(J).GE.14.0 .AND. D(J).LT.16.0) DI(J)=15 IF(D(J).GE.16.0 .AND. D(J).LT.18.0) DI(J)=17 IF(D(J).GE.18.0 .AND. D(J).LT.20.0) DI(J)=19 IF(D(J).GE.20.0 .AND. D(J).LT.22.0) DI(J)=21 IF(D(J).GE.22.0 .AND. D(J).LT.24.0) DI(J)=23 IF(D(J).GE.24.0 .AND. D(J).LT.26.0) DI(J)=25 IF(D(J).GE.26.0 .AND. D(J).LT.28.0) DI(J)=27 IF(D(J).GE.28.0 .AND. D(J).LT.30.0) DI(J)=29 IF(D(J).GE.30.0 .AND. D(J).LT.32.0) DI(J)=31 IF(D(J).GE.32.0 .AND. D(J).LT.34.0) DI(J)=33 IF(D(J).GE.34.0 .AND. D(J).LT.36.0) DI(J)=35 IF(D(J).GE.36.0 .AND. D(J).LT.38.0) DI(J)=37 IF(D(J).GE.38.0 .AND. D(J).LT.40.0) DI(J)=39 IF(D(J).GE.40.0 .AND. D(J).LT.42.0) DI(J)=41 IF(D(J).GE.42.0 .AND. D(J).LT.44.0) DI(J)=43 IF(D(J).GE.44.0 .AND. D(J).LT.46.0) DI(J)=45 IF(D(J).GE.46.0 .AND. D(J).LT.48.0) DI(J)=47 IF(D(J).GE.48.0 .AND. D(J).LT.50.0) DI(J)=49 IF(D(J).GE.50.0 .AND. D(J).LT.52.0) DI(J)=51 IF(D(J).GE.52.0 .AND. D(J).LT.54.0) DI(J)=53 IF(D(J).GE.54.0 .AND. D(J).LT.56.0) DI(J)=55 IF(D(J).GE.56.0 .AND. D(J).LT.58.0) DI(J)=57 IF(D(J).GE.58.0 .AND. D(J).LT.60.0) DI(J)=59 IF(D(J).GE.60.0 .AND. D(J).LT.62.0) DI(J)=61 IF(D(J).GE.62.0 .AND. D(J).LT.64.0) DI(J)=63 IF(D(J).GE.64.0 .AND. D(J).LT.66.0) DI(J)=65 IF(D(J).GE.66.0) DI(J)=67 IF(IPL(J).EQ.1) IIPL(J)=1 IF(IPL(J).EQ.2) IIPL(J)=2 IF(IPL(J).GE.3 .AND. IPL(J).LE.7) IIPL(J)=3 IF(IPL(J).EQ.8 .AND. IPL(J).EQ.9) IIPL(J)=2 1 CONTINUE C Tehd„„n PUUSTO.HAR (21) ja *.LEI (20) tiedostot C:\HARSI hakemistoon INQUIRE (FILE='C:\LEIMAUS\PUUSTO.HAR',EXIST=LINQ) IF (.NOT. LINQ) OPEN (21,FILE='C:\LEIMAUS\PUUSTO.HAR', + STATUS='NEW') IF (LKML.EQ.0) GOTO 9 DO 7 K=1,LKML WRITE (20,201) LEIMA(K) 7 CONTINUE 201 FORMAT (I4) 9 CONTINUE CLOSE (20) IF (LINQ) GOTO 95 IXMAX=300 IYMAX=300 IXMIN=0 IYMIN=0 WRITE (21,203) IXMIN,IYMIN,IXMAX,IYMAX 203 FORMAT (I1,1X,I1,1X,I4,1X,I4) DO 8 K=1,I WRITE (21,202) NRO(K),IIPL(K),X(K),Y(K),DI(K),PIIT(K) 202 FORMAT (I4,1X,I1,1X,I4,1X,I4,1X,I2,1X,F4.1) 8 CONTINUE CLOSE (21) 95 CONTINUE RETURN END C Subroutine aliohjelma koealan valtapituuden laskemiseksi C Laskenta ha:lla 100 paksuimman puun pituuden aritmeettisena keskiarvona C C Pituusk„yr„n H=A*(d/(d+B)) parametrit B10 ja B11 1.jakson pituuskoepuille. C koealalta (PALA) pinta-alan mukaiset AINT(100.*PALA) paksuinta puuta keski- C arvon laskentaan. C C Haetaan AINT(100.*PALA) kpl suurimmat l„pimitat D-taulukosta C I puiden lukum„„r„ koealalla SUBROUTINE VALTAP(D,I,PALA,B10,B11,HDOM,JHDOM,LEIMA,LKML) PARAMETER (MAXP=500) REAL D(MAXP),PALA,B10,B11,HDOM,JHDOM,SUMD,MEAND,APUD(MAXP) INTEGER I,J,K,L,M,KPL,LEIMA(MAXP),LKML CHARACTER APU1*1 KPL=INT(100.*PALA) DO 1 J=1,I ISO=0 DO 2 K=1,I IF(J.EQ.K) GOTO 2 IF(D(J).GT.D(K)) ISO=ISO+1 IF(ISO.GT.(I-KPL-1)) THEN L=L+1 SUMD=SUMD+D(J) MEAND=SUMD/REAL(L) GOTO 1 END IF 2 CONTINUE 1 CONTINUE HDOM=1.3 + (MEAND**2) /(B10 + B11 * MEAND)**2 C J„„v„n puuston valtapituus laskemiseksi t„ytet„„n APUD-taulukko C j„„vien puiden l„pimitoilla ja lasketaan HDOM kuten edell„, tosin C AINT(100.*PALA)*j„„vien m„„r„st„ paksuimpia puita DO 3 IA=1,I DO 4 IB=1,LKML IF(LEIMA(IB).EQ.IA) GOTO 3 4 CONTINUE M=M+1 APUD(M)=D(IA) 3 CONTINUE C Taulukossa APUD M kappaletta j„„vien puiden l„pimittoja ISO=0 SUMD=0. MEAND=0. L=0 KPL=INT(100.*PALA) DO 5 J=1,M ISO=0 DO 6 K=1,M IF(J.EQ.K) GOTO 6 IF(APUD(J).GT.APUD(K)) ISO=ISO+1 IF(ISO.GT.(M-KPL-1)) THEN L=L+1 SUMD=SUMD+APUD(J) MEAND=SUMD/REAL(L) GOTO 5 END IF 6 CONTINUE 5 CONTINUE JHDOM=1.3 + (MEAND**2)/(B10+B11*MEAND)**2 RETURN END C Marv1 Systarv. Pituusk„yrien analyyttinen laadinta C C N„slund'n pituusk„yr„ SUBROUTINE PITU (IALKU,JNRO,DI,PI,B0,B1) INTEGER I,J,K,L,IX,JX,MAXP,IALKU,JNRO REAL DI(500), PI(500), SELI(500), B0, B1, COV SUMD2=0. SUMD=0. SUMSEL=0. DO 10 I=1,JNRO IF ((PI(i)-1.3).LT.0.1) write(6,*) i, PI(i) SELI(I) = DI(I) / SQRT (PI(I)-1.3) SUMD2 = SUMD2 + DI(I)**2 SUMD = SUMD + DI(I) SUMSEL = SUMSEL + SELI(I) SUMTUL = SUMTUL + SELI(I)*DI(I) 10 CONTINUE COV = SUMTUL - SUMD*SUMSEL/REAL(JNRO) VAR = SUMD2 - SUMD**2/REAL(JNRO) B1 = COV/VAR B0 = SUMSEL/REAL(JNRO) - B1*SUMD/REAL(JNRO) RETURN END