herest.c

HMM的一种经典的训练算法,需要的快下啊

         if (NextArg()!=STRINGARG)
            HError(2319,"HERest: input transform directory expected");
         AddInXFormDir(&hset,GetStrArg());
         if (NextArg()==STRINGARG) {
            if (xfInfo.inXFormExt == NULL)
               xfInfo.inXFormExt = GetStrArg(); 
            else
               HError(2319,"HERest: only one input transform extension may be specified");
         }
	 if (NextArg() != SWITCHARG)
	   HError(2319,"HERest: cannot have -J as the last option");	  
         break;              
      case 'K':
         if (NextArg()!=STRINGARG)
            HError(2319,"HERest: output transform directory expected");
         xfInfo.outXFormDir = GetStrArg(); 
         if (NextArg()==STRINGARG)
	   xfInfo.outXFormExt = GetStrArg(); 
	 if (NextArg() != SWITCHARG)
	   HError(2319,"HERest: cannot have -K as the last option");	  
         break;              
      case 'z':
         if (NextArg() != STRINGARG)
            HError(2319,"HERest: output TMF file expected");
         xfInfo.xformTMF = GetStrArg(); break;
      default:
         HError(2319,"HERest: Unknown switch %s",s);
      }
   } 
   if (NextArg() != STRINGARG)
      HError(2319,"HERest: file name of vocabulary list expected");

   Initialise(fbInfo, &fbInfoStack, &hset, GetStrArg());
   InitUttInfo(utt, twoDataFiles);
   numUtt = 1;

   if (trace&T_TOP) 
      SetTraceFB(); /* allows HFB to do top-level tracing */

   do {
      if (NextArg()!=STRINGARG)
         HError(2319,"HERest: data file name expected");
      if (twoDataFiles && (parMode!=0)){
         if ((NumArgs() % 2) != 0)
            HError(2319,"HERest: Must be even num of training files for single pass training");
         strcpy(datafn1,GetStrArg());
         datafn = datafn1;
         
         datafn2 = GetStrArg();
      }else
         datafn = GetStrArg();
      if (parMode==0){
         src=LoadAccs(&hset, datafn,uFlags);
         ReadFloat(&src,&tmpFlt,1,ldBinary);
         totalPr += (LogDouble)tmpFlt;
         ReadInt(&src,&tmpInt,1,ldBinary);
         totalT += tmpInt;
         CloseSource( &src );
      }
      else {
         /* track speakers */	 
         if (UpdateSpkrStats(&hset,&xfInfo, datafn)) spUtt=0;
	 /* Check to see whether set-up is valid */
	 CheckUpdateSetUp();
         fbInfo->inXForm = xfInfo.inXForm;
         fbInfo->al_inXForm = xfInfo.al_inXForm;
         fbInfo->paXForm = xfInfo.paXForm;
         if ((maxSpUtt==0) || (spUtt<maxSpUtt))
            DoForwardBackward(fbInfo, utt, datafn, datafn2) ;
         numUtt += 1; spUtt++;
      }
   } while (NumArgs()>0);

   if (uFlags&UPXFORM) {/* ensure final speaker correctly handled */ 
      UpdateSpkrStats(&hset,&xfInfo, NULL); 
      if (trace&T_TOP) {
         printf("Reestimation complete - average log prob per frame = %e (%d frames)\n",
                totalPr/totalT, totalT);
      }
   } else {
      if (parMode>0  || (parMode==0 && (updateMode&UPMODE_DUMP))){
         MakeFN("HER$.acc",newDir,NULL,newFn);
         f=DumpAccs(&hset,newFn,uFlags,parMode);
         tmpFlt = (float)totalPr;
         WriteFloat(f,&tmpFlt,1,ldBinary);
         WriteInt(f,(int*)&totalT,1,ldBinary);
         fclose( f );
      }
      if (parMode <= 0) {
         if (stats) {
            StatReport(&hset);
         }
         if (updateMode&UPMODE_UPDATE)
            UpdateModels(&hset,utt->pbuf2);
      }
   }
   ResetHeap(&uttStack);
   ResetHeap(&fbInfoStack);
   ResetHeap(&hmmStack);
   Exit(0);
   return (0);          /* never reached -- make compiler happy */
}

/* -------------------------- Initialisation ----------------------- */

void Initialise(FBInfo *fbInfo, MemHeap *x, HMMSet *hset, char *hmmListFn)
{   
   HSetKind hsKind;
   int L,P,S,vSize,maxM; 

   /* Load HMMs and init HMMSet related global variables */
   if(MakeHMMSet( hset, hmmListFn )<SUCCESS)
      HError(2321,"Initialise: MakeHMMSet failed");
   if(LoadHMMSet( hset,hmmDir,hmmExt)<SUCCESS)
      HError(2321,"Initialise: LoadHMMSet failed");
   AttachAccs(hset, &accStack, uFlags);
   ZeroAccs(hset, uFlags);
   P = hset->numPhyHMM;
   L = hset->numLogHMM;
   vSize = hset->vecSize;
   S = hset->swidth[0];
   maxM = MaxMixInSet(hset);

   hsKind = hset->hsKind;
   if (hsKind==DISCRETEHS)
     uFlags = (UPDSet) (uFlags & (~(UPMEANS|UPVARS|UPXFORM)));

   if (parMode != 0) {
      ConvDiagC(hset,TRUE);
   }
   if (trace&T_TOP) {
      if (uFlags&UPMAP)  printf("HERest  MAP Updating: ");
      else printf("HERest  ML Updating: ");
      if (uFlags&UPTRANS) printf("Transitions "); 
      if (uFlags&UPMEANS) printf("Means "); 
      if (uFlags&UPVARS)  printf("Variances "); 
      if (uFlags&UPXFORM)  printf("XForms "); 
      if (uFlags&UPMIXES && maxM>1)  printf("MixWeights "); 
      printf("\n\n ");
    
      if (parMode>=0) printf("Parallel-Mode[%d] ",parMode);

      printf("System is ");
      switch (hsKind){
      case PLAINHS:  printf("PLAIN\n");  break;
      case SHAREDHS: printf("SHARED\n"); break;
      case TIEDHS:   printf("TIED\n"); break;
      case DISCRETEHS: printf("DISCRETE\n"); break;
      }

      printf("%d Logical/%d Physical Models Loaded, VecSize=%d\n",L,P,vSize);
      if (hset->numFiles>0)
         printf("%d MMF input files\n",hset->numFiles);
      if (mmfFn != NULL)
         printf("Output to MMF file:  %s\n",mmfFn); 
      fflush(stdout);
   }
   SetVFloor( hset, vFloor, minVar);
   totalPr = 0.0;

   if (xfInfo.inSpkrPat == NULL) xfInfo.inSpkrPat = xfInfo.outSpkrPat; 
   if (xfInfo.paSpkrPat == NULL) xfInfo.paSpkrPat = xfInfo.outSpkrPat; 
   if (uFlags&UPXFORM) {
      if ((hsKind != PLAINHS) && (hsKind != SHAREDHS))
         HError(999,"Can only estimated transforms with PLAINHS and SHAREDHS!");
      if (uFlags != UPXFORM)
         HError(999,"Can only update linear transforms OR model parameters!");
      xfInfo.useOutXForm = TRUE;
      /* This initialises things - temporary hack - THINK!! */
      CreateAdaptXForm(hset, "tmp");
   } 

   
   /* initialise and  pass information to the forward backward library */
   InitialiseForBack(fbInfo, x, hset, uFlags, pruneInit, pruneInc,
                     pruneLim, minFrwdP);

   if (parMode != 0) {
      ConvLogWt(hset);
   }
   /* 2-model reestimation */
   if (al_hmmUsed){
       if (trace&T_TOP)
           printf("2-model re-estimation enabled\n");
       /* load alignment HMM set */
       CreateHMMSet(&al_hset,&hmmStack,TRUE);
       xfInfo.al_hset = &al_hset;
       if (xfInfo.alXFormExt == NULL) xfInfo.alXFormExt = xfInfo.inXFormExt;
       /* load multiple MMFs */
       if (strlen(al_hmmMMF) > 0 ) {
           char *p,*q;
           Boolean eos;
           p=q=al_hmmMMF;
           for(;;) {
               eos = (*p=='\0');
               if ( ( isspace((int) *p) || *p == '\0' ) && (q!=p) ) {
                   *p='\0';
                   if (trace&T_TOP) { 
                       printf("Loading alignment HMM set %s\n",q);
                   }
                   AddMMF(&al_hset,q);
                   if (eos)
                       break;
                   q=p+1;
               }
               p++;
           }
       }
       if (strlen(al_hmmLst) > 0 ) 
           MakeHMMSet(&al_hset, al_hmmLst );
       else /* use same hmmList */
           MakeHMMSet(&al_hset, hmmListFn );
       if (strlen(al_hmmDir) > 0 )
           LoadHMMSet(&al_hset,al_hmmDir,al_hmmExt);
       else
           LoadHMMSet(&al_hset,NULL,NULL);

       /* switch model set */
       UseAlignHMMSet(fbInfo,x,&al_hset);
       if (parMode != 0) {
	  ConvDiagC(&al_hset,TRUE);
	  ConvLogWt(&al_hset);
       }

       /* and echo status */
       if (trace&T_TOP) { 
           if (strlen(al_hmmDir) > 0 )
               printf(" HMM Dir %s",al_hmmDir);
           if (strlen(al_hmmExt) > 0 )
               printf(" Ext %s",al_hmmExt);
           printf("\n");
           if (strlen(al_hmmLst) > 0 )
               printf("HMM List %s\n",al_hmmLst);
           printf(" %d Logical/%d Physical Models Loaded, VecSize=%d\n",
                  al_hset.numLogHMM,al_hset.numPhyHMM,al_hset.vecSize);
       }
   }
}

/* ------------------- Statistics Reporting  -------------------- */

/* PrintStats: for given hmm */
void PrintStats(HMMSet *hset,FILE *f, int n, HLink hmm, int numEgs)
{
   WtAcc *wa;
   char buf[MAXSTRLEN];
   StateInfo *si;
   int i,N;
    
   N = hmm->numStates;
   ReWriteString(HMMPhysName(hset,hmm),buf,DBL_QUOTE);
   fprintf(f,"%4d %14s %4d ",n,buf,numEgs);
   for (i=2;i<N;i++) {
      si = hmm->svec[i].info;
      wa = (WtAcc *)((si->pdf+1)->hook);
      fprintf(f," %10f",wa->occ);
   }
   fprintf(f,"\n");
}

/* StatReport: print statistics report */
void StatReport(HMMSet *hset)
{
   HMMScanState hss;
   HLink hmm;
   FILE *f;
   int px;

   if ((f = fopen(statFN,"w")) == NULL){
      HError(2311,"StatReport: Unable to open stats file %s",statFN);
      return;
   }
   NewHMMScan(hset,&hss);
   px=1;
   do {
      hmm = hss.hmm;
      PrintStats(hset,f,px,hmm,(int)hmm->hook);
      px++;
   } while (GoNextHMM(&hss));
   EndHMMScan(&hss);
   fclose(f);
}

/* -------------------- Top Level of F-B Updating ---------------- */


/* Load data and call FBFile: apply forward-backward to given utterance */
void DoForwardBackward(FBInfo *fbInfo, UttInfo *utt, char * datafn, char * datafn2)
{
   utt->twoDataFiles = twoDataFiles ;
   utt->S = fbInfo->al_hset->swidth[0];

   /* Load the labels */
   LoadLabs(utt, lff, datafn, labDir, labExt);
   /* Load the data */
   LoadData(fbInfo->al_hset, utt, dff, datafn, datafn2);

   if (firstTime) {
      InitUttObservations(utt, fbInfo->al_hset, datafn, fbInfo->maxMixInS);
      firstTime = FALSE;
   }
  
   /* fill the alpha beta and otprobs (held in fbInfo) */
   if (FBFile(fbInfo, utt, datafn)) {
      /* update totals */
      totalT += utt->T ;
      totalPr += utt->pr ;
   }
}

/* --------------------------- Model Update --------------------- */

static int nFloorVar = 0;     /* # of floored variance comps */
static int nFloorVarMix = 0;  /* # of mix comps with floored vars */

/* UpdateTrans: use acc values to calc new estimate for transP */
void UpdateTrans(HMMSet *hset, int px, HLink hmm)
{
   int i,j,N;
   float x,occi;
   TrAcc *ta;
   
   ta = (TrAcc *) GetHook(hmm->transP);
   if (ta==NULL) return;   /* already done */
   N = hmm->numStates;
   for (i=1;i<N;i++) {
      occi = ta->occ[i];
      if (occi > 0.0) 
         for (j=2;j<=N;j++) {
            x = ta->tran[i][j]/occi;
            hmm->transP[i][j] = (x>MINLARG)?log(x):LZERO;
         }
      else
         HError(-2326,"UpdateTrans: Model %d[%s]: no transitions out of state %d",
                px,HMMPhysName(hset,hmm),i);
   }
   SetHook(hmm->transP,NULL);
}

/* FloorMixes: apply floor to given mix set */
void FloorMixes(HMMSet *hset, MixtureElem *mixes, int M, float floor)
{
   float sum,fsum,scale;
   MixtureElem *me;
   int m;
   
   if (hset->logWt == TRUE) HError(999,"FloorMixes requires linear weights");
   sum = fsum = 0.0;
   for (m=1,me=mixes; m<=M; m++,me++) {
      if (MixWeight(hset,me->weight)>floor)
         sum += me->weight;
      else {
         fsum += floor; me->weight = floor;
      }
   }
   if (fsum>1.0) HError(2327,"FloorMixes: Floor sum too large");
   if (fsum == 0.0) return;
   if (sum == 0.0) HError(2328,"FloorMixes: No mixture weights above floor");
   scale = (1.0-fsum)/sum;
   for (m=1,me=mixes; m<=M; m++,me++)
      if (me->weight>floor) me->weight *= scale;
}

/* FloorTMMixes: apply floor to given tied mix set */
void FloorTMMixes(Vector mixes, int M, float floor)
{
   float sum,fsum,scale,fltWt;
   int m;
   
   sum = fsum = 0.0;
   for (m=1; m<=M; m++) {
      fltWt = mixes[m];
      if (fltWt>floor)
         sum += fltWt;
      else {
         fsum += floor;
         mixes[m] = floor;
      }
   }
   if (fsum>1.0) HError(2327,"FloorTMMixes: Floor sum too large");
   if (fsum == 0.0) return;
   if (sum == 0.0) HError(2328,"FloorTMMixes: No mixture weights above floor");
   scale = (1.0-fsum)/sum;
   for (m=1; m<=M; m++){
      fltWt = mixes[m];
      if (fltWt>floor)
         mixes[m] = fltWt*scale;
   }
}

/* FloorDProbs: apply floor to given discrete prob set */
void FloorDProbs(ShortVec mixes, int M, float floor)
{
   float sum,fsum,scale,fltWt;
   int m;
   
   sum = fsum = 0.0;
   for (m=1; m<=M; m++) {
      fltWt = Short2DProb(mixes[m]);
      if (fltWt>floor)
         sum += fltWt;
      else {
         fsum += floor;
         mixes[m] = DProb2Short(floor);
      }
   }
   if (fsum>1.0) HError(2327,"FloorDProbs: Floor sum too large");
   if (fsum == 0.0) return;
   if (sum == 0.0) HError(2328,"FloorDProbs: No probabilities above floor");
   scale = (1.0-fsum)/sum;
   for (m=1; m<=M; m++){
      fltWt = Short2DProb(mixes[m]);
      if (fltWt>floor)
         mixes[m] = DProb2Short(fltWt*scale);
   }
}

/* UpdateWeights: use acc values to calc new estimate of mix weights */
void UpdateWeights(HMMSet *hset, int px, HLink hmm)
{
   int i,s,m,M=0,N,S;
   float x,occi;
   WtAcc *wa;
   StateElem *se;
   StreamElem *ste;
   MixtureElem *me;
   HSetKind hsKind;

   N = hmm->numStates;
   se = hmm->svec+2;
   hsKind = hset->hsKind;
   S = hset->swidth[0];
   for (i=2; i<N; i++,se++){
      ste = se->info->pdf+1;
      for (s=1;s<=S; s++,ste++){
         wa = (WtAcc *)ste->hook;
         switch (hsKind){
         case TIEDHS:
            M=hset->tmRecs[s].nMix;
            break;
         case DISCRETEHS:
         case PLAINHS:
         case SHAREDHS:
            M=ste->nMix;
            break;
         }
         if (wa != NULL) {
            occi = wa->occ;
            if (occi>0) {
               for (m=1; m<=M; m++){
                  x = wa->c[m]/occi;
                  if (x>1.0){
                     if (x>1.001)