I am struggling to replicate the following fortran 77 subroutine into C# method. The main problem is that I am trying to get rid of the fortran's goto statements.
The fortran to C# conversion looks fine, still the goto statements remained. Is there a way all of them can be replaced with conditionals or some other ways? Thank you.
Here is the fortran subroutine:
C PROGRAM TEMP
      subroutine TEMP (acl,adu,aeff,c,cair,cb,cbare, 
     + cclo,count1,csum,di,ed,emcl,emsk,enbal, 
     + enbal2,ere,erel,esw,eswdif,eswphy,eswpot, 
     + evap,facl,fcl,fec,feff,food,h,hc,he,ht,htcl,icl,j, 
     + mbody,p,po,r1,r2,rbare,rcl, 
     + rclo,rclo2,rdcl,rdsk,rob,rsum,sex,sigm,sw,swf,swm, 
     + ta,tbody,tcl,tcore,tmrt,tsk,v,vb,vb1,vb2, 
     + vpa,vpts,wetsk,wd,wr,ws,wsum,xx) 
      real acl,adu,aeff,c(0:10),cair,cb,cbare, 
     + cclo,csum,di,ed,emcl,emsk,enbal, 
     + enbal2,ere,erel,esw,eswdif,eswphy,eswpot, 
     + evap,facl,fcl,fec,feff,food,h,hc,he,ht,htcl,icl, 
     + mbody,p,po,r1,r2,rbare,rcl, 
     + rclo,rclo2,rdcl,rdsk,rob,rsum,sigm,sw,swf,swm, 
     + ta,tbody,tcl,tcore(1:7),tmrt,tsk,v,vb,vb1,vb2, 
     + vpa,vpts,wetsk,wd,wr,ws,wsum,xx 
      integer count1,count3,j,sex 
      wetsk = 0. 
      adu = 0.203 * mbody ** 0.425 * ht ** 0.725 
      hc = 2.67 + ( 6.5 * v ** 0.67) 
      hc = hc * (p /po) ** 0.55 
      feff = 0.725 
C 
      facl = (- 2.36 + 173.51 * icl - 100.76 * icl * icl + 19.28 
     + * (icl ** 3.)) / 100. 
C 
      if (facl .gt.1.) facl = 1. 
      rcl = (icl/6.45)/facl 
      if (icl.ge.2.) y = 1. 
      if ((icl .gt. 0.6) .and. (icl .lt. 2.)) y = (ht - 0.2) / ht 
      if ((icl .le. 0.6) .and. (icl .gt. 0.3)) y = 0.5 
      if ((icl .le. 0.3) .and. (icl .gt. 0.)) y = 0.1 
      r2 = adu * (fcl - 1. + facl) / (2. * 3.14 * ht * y) 
      r1 = facl * adu / (2. * 3.14 * ht * y) 
      di = r2 - r1 
C TEMPERATURE
      do 90 j = 1,7 
      tsk = 34. 
      count1 = 0 
      tcl = (ta + tmrt + tsk) / 3. 
      count3 = 1 
      enbal2 = 0. 
   20 acl = adu * facl + adu * (fcl - 1.) 
      rclo2 = emcl*sigm *((tcl+273.2)** 4.-(tmrt+273.2)** 4.)*feff 
      htcl = 6.28 * ht * y * di / (rcl * alog(r2/r1) * acl) 
      tsk = 1. / htcl * (hc * (tcl - ta) + rclo2) + tcl 
C RADIATION 
      aeff = adu * feff 
      rbare = aeff * (1.-facl) * emsk * sigm * 
     + ((tmrt + 273.2) ** 4. - (tsk + 273.2) ** 4.) 
      rclo = feff * acl * emcl * sigm * 
     + ((tmrt + 273.2) ** 4. - (tcl + 273.2) ** 4.) 
      rsum = rbare + rclo 
C CONVECT
      cbare = hc * (ta - tsk) * adu * (1. - facl) 
      cclo = hc * (ta - tcl ) * acl 
      csum = cbare + cclo 
C CORE
      c(0) = h + ere 
      c(1) = adu * rob * cb 
      c(2) = 18. - 0.5 * tsk 
      c(3) = 5.28 * adu * c(2) 
      c(4) = 0.0208 * c(1) 
      c(5) = 0.76075 * c(1) 
      c(6) = c(3) - c(5) - tsk * c(4) 
      c(7) = - c(0) * c(2) - tsk * c(3) + tsk * c(5) 
      c(8) = c(6) * c(6) - 4. * c(4) * c(7) 
      c(9) = 5.28 * adu - c(5) - c(4) * tsk 
      c(10) = c(9) * c(9) - 4. * c(4) * 
     + (c(5) * tsk - c(0) - 5.28 * adu * tsk) 
C 
      if (tsk.eq.36.) tsk=36.01 
      tcore(7) = c(0) / (5.28 * adu + c(1) * 6.3 / 3600.) + tsk 
      tcore(3) = c(0) / (5.28 * adu + (c(1) * 6.3 / 3600.) / 
     + (1 + 0.5 * (34. -tsk))) + tsk 
      if (c(10) .lt. 0.) goto 22 
      tcore(6) = (- c(9) - c(10) ** 0.5) / (2. * c(4)) 
      tcore(1) = (- c(9) + c(10) ** 0.5) / (2. * c(4)) 
   22 if (c(8) .lt. 0.) goto 24 
      tcore(2) = (- c(6) + abs(c(8)) ** 0.5) / (2. * c(4)) 
      tcore(5) = (- c(6) - abs(c(8)) ** 0.5) / (2. * c(4)) 
   24 tcore(4) = c(0) / (5.28 * adu + c(1) * 1. / 40.) + tsk 
C TRANSPARENCE
      tbody = 0.1 * tsk + 0.9 * tcore (j) 
      swm = 304.94 * (tbody - 36.6) * adu / 3600000. 
      vpts = 6.11 * 10. ** (7.45 * tsk / (235. + tsk)) 
      if (tbody .le. 36.6) swm = 0. 
      swf = 0.7 * swm 
      if(sex .eq. 1) sw = swm 
      if(sex .eq. 2) sw = swf 
      eswphy = - sw * evap 
      he = 0.633 * hc / (p * cair) 
      fec = 1. / (1. + 0.92 * hc * rcl) 
      eswpot = he * (vpa - vpts) * adu * evap * fec 
      wetsk = eswphy / eswpot 
      if (wetsk .gt. 1.) wetsk = 1. 
      eswdif = eswphy - eswpot 
      if (eswdif .le. 0.) esw = eswpot 
      if (eswdif .gt. 0.) esw = eswphy 
      if (esw .gt. 0.) esw = 0. 
C DIFFERENCE
      rdsk = 0.79 * 10. ** 7. 
      rdcl = 0. 
      ed = evap / (rdsk + rdcl) * adu * (1 - wetsk) * (vpa-vpts) 
C VB 
      vb1 = 34. - tsk 
      vb2 = tcore(j) - 36.6 
      if (vb2 .lt.0.) vb2 = 0. 
      if (vb1 .lt.0.) vb1 = 0. 
      vb = (6.3 + 75. * (vb2)) / (1. + 0.5 * vb1) 
C BALANCE
      enbal = h + ed + ere + esw + csum + rsum + food 
C COVER
      if (count1 .eq.0) xx = 1. 
      if (count1 .eq.1) xx = 0.1 
      if (count1 .eq.2) xx = 0.01 
      if (count1 .eq.3) xx = 0.001 
      if (enbal .gt. 0.) tcl = tcl + xx 
      if (enbal .lt. 0.) tcl = tcl - xx 
      if ((enbal .le. 0.) .and. (enbal2 .gt. 0.)) goto 30 
      if ((enbal .ge. 0.) .and. (enbal2 .lt. 0.)) goto 30 
      enbal2 = enbal 
      count3 = count3 + 1 
C 
      if (count3 .gt. 200) goto 30 
      goto 20 
   30 if ((count1 .eq.0.).or.(count1.eq.1.).or.(count1.eq.2.)) then 
         count1 = count1 + 1. 
         enbal2 = 0. 
         goto 20 
      end if 
C 
      if (count1 .eq. 3.) then 
C 
         if ((j .eq. 2) .or. (j .eq. 5)) goto 40 
         if ((j .eq. 6) .or. (j .eq. 1)) goto 50 
         if (j .eq. 3) goto 60 
         if (j .eq. 7) goto 70 
         if (j .eq. 4) goto 80 
      end if 
   40 if (c(8) .lt. 0.) goto 90 
      if ((tcore(j) .ge. 36.6) .and. (tsk .le. 34.050)) goto 80 
      goto 90 
   50 if (c(10) .lt. 0. ) goto 90 
      if ((tcore(j) .ge. 36.6) .and. (tsk .gt. 33.850)) goto 80 
      goto 90 
   60 if ((tcore(j) .lt. 36.6) .and. (tsk .le. 34.000)) goto 80 
      goto 90 
   70 if ((tcore(j) .lt. 36.6) .and. (tsk .gt. 34.000)) goto 80 
      goto 90 
   80 if ((j .ne. 4) .and. (vb .ge. 91.)) goto 90 
      if ((j. eq. 4) .and. (vb .lt. 89.)) goto 90 
      if (vb .gt. 90.) vb = 90. 
C LOSSES
      ws = sw * 3600. * 1000. 
      if (ws .gt.2000.) ws = 2000. 
      wd = ed / evap * 3600. * (-1000.) 
      wr = erel / evap * 3600. * (-1000.) 
      wsum = ws + wr + wd 
      goto 100 
   90 continue 
  100 return 
      end 
And here is the C# attempt to replicate the upper fortran subroutine:
using System;
public static class GlobalMembers_TEMP
{
//C PROGRAM TEMP
// --------------------------------------------
    public static void TEMP(acl, adu, aeff, c, cair, cb, cbare, cclo, count1, csum, di, ed , emcl , emsk , enbal , enbal2 , ere , erel , esw , eswdif , eswphy , eswpot , evap , facl , fcl , fec , feff , food , h , hc , he , ht , htcl , icl , j , mbody , p , po , r1 , r2 , rbare , rcl , rclo , rclo2 , rdcl , rdsk , rob , rsum , sex , sigm , sw , swf , swm , ta , tbody , tcl, tcore , tmrt , tsk , v , vb , vb1 , vb2 , vpa , vpts , wetsk , wd , wr , ws , wsum , xx)
    {
        float acl;
        float adu;
        float aeff;
        float[] c = new float[11];  // ???
        float cair;
        float cb;
        float cbare;
        float cclo;
        float csum;
        float di;
        float ed;
        float emcl;
        float emsk;
        float enbal;
        float enbal2;
        float ere;
        float erel;
        float esw;
        float eswdif;
        float eswphy;
        float eswpot;
        float evap;
        float facl;
        float fcl;
        float fec;
        float feff;
        float food;
        float h;
        float hc;
        float he;
        float ht;
        float htcl;
        float icl;
        float mbody;
        float p;
        float po;
        float r1;
        float r2;
        float rbare;
        float rcl;
        float rclo;
        float rclo2;
        float rdcl;
        float rdsk;
        float rob;
        float rsum;
        float sigm;
        float sw;
        float swf;
        float swm;
        float ta;
        float tbody;
        float tcl;
        float[] tcore = new float[7];   // ???
        float tmrt;
        float tsk;
        float v;
        float vb;
        float vb1;
        float vb2;
        float vpa;
        float vpts;
        float wetsk;
        float wd;
        float wr;
        float ws;
        float wsum;
        float xx;
        int count1;
        int count3;
        int j;
        int sex;
        wetsk = 0.0;
        adu = 0.203 * mbody * *0.425 * ht * *0.725;
        hc = 2.67 + (6.5 * v * *0.67);
        hc = hc * (p / po) * *0.55;
        feff = 0.725;
        //C rcl = icl / 6.45
        facl = (-2.36 + 173.51 * icl - 100.76 * icl * icl + 19.28 * (icl * *3.0)) / 100.0;
        //C
        if (facl > 1.0)
        {
            facl = 1.0;
        }
        rcl = (icl / 6.45) / facl;
        if (icl >= 2.0)
        {
            y = 1.0;
        }
        if ((icl > 0.6)) && ((icl < 2.0)) 
        {
            y = (ht - 0.2) / ht;
        }
        if ((icl <= 0.6)) && ((icl > 0.3)) 
        {
            y = 0.5;
        }
        if ((icl <= 0.3)) && ((icl > 0.0)) 
        {
            y = 0.1;
        }
        r2 = adu * (fcl - 1.0 + facl) / (2.0 * 3.14 * ht * y);
        r1 = facl * adu / (2.0 * 3.14 * ht * y);
        di = r2 - r1;
        //C TEMPERATURE
        for (90 j = 1, 7)
        {
            tsk = 34.0;
            count1 = 0;
            tcl = (ta + tmrt + tsk) / 3.0;
            count3 = 1;
            enbal2 = 0.0;
    g20:
            acl = adu * facl + adu * (fcl - 1.0);
            rclo2 = emcl * sigm * ((tcl + 273.2) * *4.0 - (tmrt + 273.2) * *4.0) * feff;
            htcl = 6.28 * ht * y * di / (rcl * alog(r2 / r1) * acl);
            tsk = 1.0 / htcl * (hc * (tcl - ta) + rclo2) + tcl;
            //C RADIATION
            aeff = adu * feff;
            rbare = aeff * (1.0 - facl) * emsk * sigm * ((tmrt + 273.2) * *4.0 - (tsk + 273.2) * *4.0);
            rclo = feff * acl * emcl * sigm * ((tmrt + 273.2) * *4.0 - (tcl + 273.2) * *4.0);
            rsum = rbare + rclo;
            //C CONVECT
            cbare = hc * (ta - tsk) * adu * (1.0 - facl);
            cclo = hc * (ta - tcl) * acl;
            csum = cbare + cclo;
            //C CORE
            c[0] = h + ere;
            c[1] = adu * rob * cb;
            c[2] = 18.0 - 0.5 * tsk;
            c[3] = 5.28 * adu * c[2];
            c[4] = 0.0208 * c[1];
            c[5] = 0.76075 * c[1];
            c[6] = c[3] - c[5] - tsk * c[4];
            c[7] = -c[0] * c[2] - tsk * c[3] + tsk * c[5];
            c[8] = c[6] * c[6] - 4.0 * c[4] * c[7];
            c[9] = 5.28 * adu - c[5] - c[4] * tsk;
            c[10] = c[9] * c[9] - 4.0 * c[4] * (c[5] * tsk - c[0] - 5.28 * adu * tsk);
            //C
            if (tsk == 36.0)
            {
                tsk = 36.01;
            }
            tcore[7] = c[0] / (5.28 * adu + c[1] * 6.3 / 3600.0) + tsk;
            tcore[3] = c[0] / (5.28 * adu + (c[1] * 6.3 / 3600.0) / (1 + 0.5 * (34.0 - tsk))) + tsk;
            if (c[10] < 0.0)
            {
                goto g22;
            }
            tcore[6] = (-c[9] - c[10] * *0.5) / (2.0 * c[4]);
            tcore[1] = (-c[9] + c[10] * *0.5) / (2.0 * c[4]);
    g22:
            if (c[8] < 0.0)
            {
                goto g24;
            }
            tcore[2] = (-c[6] + Math.Abs(c[8]) * *0.5) / (2.0 * c[4]);
            tcore[5] = (-c[6] - Math.Abs(c[8]) * *0.5) / (2.0 * c[4]);
    g24:
            tcore[4] = c[0] / (5.28 * adu + c[1] * 1.0 / 40.0) + tsk;
            //C TRANSPARENCE
            tbody = 0.1 * tsk + 0.9 * tcore[j];
            swm = 304.94 * (tbody - 36.6) * adu / 3600000.0;
            vpts = 6.11 * 10.0 * *(7.45 * tsk / (235.0 + tsk));
            if (tbody <= 36.6)
            {
                swm = 0.0;
            }
            swf = 0.7 * swm;
            if (sex == 1)
            {
                sw = swm;
            }
            if (sex == 2)
            {
                sw = swf;
            }
            eswphy = -sw * evap;
            he = 0.633 * hc / (p * cair);
            fec = 1.0 / (1.0 + 0.92 * hc * rcl);
            eswpot = he * (vpa - vpts) * adu * evap * fec;
            wetsk = eswphy / eswpot;
            if (wetsk > 1.0)
            {
                wetsk = 1.0;
            }
            eswdif = eswphy - eswpot;
            if (eswdif <= 0.0)
            {
                esw = eswpot;
            }
            if (eswdif > 0.0)
            {
                esw = eswphy;
            }
            if (esw > 0.0)
            {
                esw = 0.0;
            }
            //C DIFFERENCE
            rdsk = 0.79 * 10.0 * *7.0;
            rdcl = 0.0;
            ed = evap / (rdsk + rdcl) * adu * (1 - wetsk) * (vpa - vpts);
            //C VB
            vb1 = 34.0 - tsk;
            vb2 = tcore[j] - 36.6;
            if (vb2 < 0.0)
            {
                vb2 = 0.0;
            }
            if (vb1 < 0.0)
            {
                vb1 = 0.0;
            }
            vb = (6.3 + 75.0 * (vb2)) / (1.0 + 0.5 * vb1);
            //C BALANCE
            enbal = h + ed + ere + esw + csum + rsum + food;
            //C COVER
            if (count1 == 0)
            {
                xx = 1.0;
            }
            if (count1 == 1)
            {
                xx = 0.1;
            }
            if (count1 == 2)
            {
                xx = 0.01;
            }
            if (count1 == 3)
            {
                xx = 0.001;
            }
            if (enbal > 0.0)
            {
                tcl = tcl + xx;
            }
            if (enbal < 0.0)
            {
                tcl = tcl - xx;
            }
            if ((enbal <= 0.0)) && ((enbal2 > 0.0)) 
            {
                goto g30;
            }
            if ((enbal >= 0.0)) && ((enbal2 < 0.0)) 
            {
                goto g30;
            }
            enbal2 = enbal;
            count3 = count3 + 1;
            //C
            if (count3 > 200)
            {
                goto g30;
            }
            goto g20;
    g30:
            if ((count1 == 0.0))
            {
                || ((count1 == 1.0)) || ((count1 == 2.0))
                {
                count1 = count1 + 1.0;
                enbal2 = 0.0;
                goto g20;
                }
            }
            //C
            if (count1 == 3.0)
            {
                //C
                if ((j == 2))
                {
                    || ((j == 5)) goto g40;
                }
                if ((j == 6))
                {
                    || ((j == 1)) goto g50;
                }
                if (j == 3)
                {
                    goto g60;
                }
                if (j == 7)
                {
                    goto g70;
                }
                if (j == 4)
                {
                    goto g80;
                }
            }
    g40:
            if (c[8] < 0.0)
            {
                goto g90;
            }
            if ((tcore[j] >= 36.6)) && ((tsk <= 34.050))
            {
                goto g80;
            }
            goto g90;
    g50:
            if (c[10] < 0.0)
            {
                goto g90;
            }
            if ((tcore[j] >= 36.6)) && ((tsk > 33.850))
            {
                goto g80;
            }
            goto g90;
    g60:
            if ((tcore[j] < 36.6)) && ((tsk <= 34.000))
            {
                goto g80;
            }
            goto g90;
    g70:
            if ((tcore[j] < 36.6)) && ((tsk > 34.000))
            {
                goto g80;   
            }
            goto g90;
    g80:
            if ((j != 4)) && ((vb >= 91.0)) 
            {
                goto g90;
            }
            if ((j.eq.4)) && ((vb < 89.0))
            {
                 goto g90;
            }
            if (vb > 90.0)
            {
                vb = 90.0;
            }
            //C LOSSES
            ws = sw * 3600.0 * 1000.0;
            if (ws > 2000.0)
            {
                ws = 2000.0;
            }
            wd = ed / evap * 3600.0 * (-1000.0);
            wr = erel / evap * 3600.0 * (-1000.0);
            wsum = ws + wr + wd;
            goto g100;
        }
    g100:
        return;
    }
}
 
     
     
     
    