/*-*-C-*-*/ program Bic_Lum /* **MEMBER**=SLCLIBS:xxxSHRLIB **CMS**=xxxUTIL ====================================================================== Abs: This sequence reads the processed bunch currents and counts filled buckets for each ring plus the overlap of filled buckets for luminosity calculation purposes Name: Bic_Lum.ST Proto: none Auth: 12-Apr-1999, Ron Chestnut Rev: DD-MMM-YYYY, Reviewer's Name (.NE. Author's Name) ----------------------------------------------------------------------- Mod: 7-May-2001 Add checks for BIT-3 memory; publish errors 12-Oct-2000 Add sq_i_l, _h, and _hl for sum over buckets of current squared 11-Nov-1999 Add LAST_3_SHIFTS for Luminosity display 26-Oct-1999 Add minute_of_day for daylight savings time correctness Move from LUMVAL to LUMCOR for corrected luminosity 2-Sep-1999 Add check for BIC downtime; zero intervening data 6-Aug-1999 Add shift totals 7-July-1999 Get SPLUMVAL caluclated every second 21-Aug-2000 Add Luminosity Lifetime Calculations =======================================================================*/ #include "bic_common.h" %%#include "math.h" %%#include "epicsPrint.h" /*epicsPrintf */ %%#include "taskLib.h" /* TaskDelay */ %%#include "tsDefs.h" /* Time functions and structures */ %%#include "string.h" /* Time functions and structures */ %%#include "stdio.h" /* Time functions and structures */ #define PUBLISH_LER(message) \ sprintf(state_msg_ler,message);pvPut(state_msg_ler); \ epicsPrintf("LER: %s\n",state_msg_ler); #define PUBLISH_HER(message) \ sprintf(state_msg_her,message);pvPut(state_msg_her); \ epicsPrintf("HER: %s\n",state_msg_her); #define N_MINUTES 1440 #define HERI_LOC 0xf002a000 #define LERI_LOC 0xf002c000 #define LUM_LOC 0xf002e000 #define INT_CONST 0.00006 /* Map to all the relevant Database items */ string state_msg_ler; assign state_msg_ler to "LB60:STATE:MSG1"; string state_msg_her; assign state_msg_her to "HB60:STATE:MSG1"; float day_heri[N_MINUTES]; assign day_heri to "PB60:DAILY:HERI"; float day_leri[N_MINUTES]; assign day_leri to "PB60:DAILY:LERI"; float day_lum[N_MINUTES]; assign day_lum to "PB60:DAILY:LUM"; float day_tmark[N_MINUTES]; assign day_tmark to "PB60:DAILY:TMARK"; float day_tmod24[N_MINUTES]; assign day_tmod24 to "PB60:DAILY:TMOD24"; float day_runavg[N_MINUTES]; assign day_runavg to "PB60:DAILY:RUNAVG"; float runavg_temp; long i_min; assign i_min to "PB60:DAILY:MINUTE"; long i_min_old; float bc_l[N_BUCKET]; assign bc_l to "LB60:BNCHCURR:GOOD"; float bc_h[N_BUCKET]; assign bc_h to "HB60:BNCHCURR:GOOD"; float bc_o[N_BUCKET]; assign bc_o to "PB60:BNCHCURR:GOOD"; float goal_l[N_BUCKET]; assign goal_l to "LB60:FILL:GOAL"; float goal_h[N_BUCKET]; assign goal_h to "HB60:FILL:GOAL"; float dcct_l; assign dcct_l to "LB60:DCCT:SUMY"; monitor dcct_l; float dcct_h; assign dcct_h to "HB60:DCCT:SUMY"; monitor dcct_h; float lumval; assign lumval to "PB60:LUMCOR"; monitor lumval; float lummed; assign lummed to "PB60:LUM:MEDFVAL"; monitor lummed; float splumval; assign splumval to "PB60:SPLUMVAL"; float splumval_inv; float lumcoll; assign lumcoll to "PB60:LUMCOR:COLL"; monitor lumcoll; float lumdaq; assign lumdaq to "PB60:LUMCOR:DAQ"; monitor lumdaq; float dcct_ler_ped; assign dcct_ler_ped to "LB60:DCCT:SUMY.A"; float dcct_her_ped; assign dcct_her_ped to "HB60:DCCT:SUMY.A"; monitor dcct_ler_ped; monitor dcct_her_ped; float ped_ler_old; float ped_her_old; float buck_lum; assign buck_lum to "PB60:LUM:BUCK"; float sq_i_l; assign sq_i_l to "LB60:SUM:ISQUARE"; float sq_i_h; assign sq_i_h to "HB60:SUM:ISQUARE"; float sq_i_hl; assign sq_i_hl to "PB60:SUM:ISQUARE"; long num_filled_l; assign num_filled_l to "LB60:BNCHCURR:FILL:NUMBUCK"; long num_filled_h; assign num_filled_h to "HB60:BNCHCURR:FILL:NUMBUCK"; long num_filled_o; assign num_filled_o to "PB60:FILL:NUMBUCK"; long num_saved_o; assign num_saved_o to "PB60:SAVE:NUMBUCK"; long num_saved_h; assign num_saved_h to "HB60:SAVE:NUMBUCK"; long num_saved_l; assign num_saved_l to "LB60:SAVE:NUMBUCK"; monitor num_saved_o; monitor num_saved_h; monitor num_saved_l; short current_shift; assign current_shift to "PB60:SHIFT"; short old_shift; float current_tot_lum; assign current_tot_lum to "PB60:SHIFT:TOT"; monitor current_tot_lum; float current_col_lum; assign current_col_lum to "PB60:SHIFT:COL"; monitor current_col_lum; float current_daq_lum; assign current_daq_lum to "PB60:SHIFT:DAQ"; monitor current_daq_lum; float best_tot_lum[3]; assign best_tot_lum to {"PB60:BTOT:OWL","PB60:BTOT:DAY","PB60:BTOT:SWING"}; float best_col_lum[3]; assign best_col_lum to {"PB60:BCOL:OWL","PB60:BCOL:DAY","PB60:BCOL:SWING"}; float best_daq_lum[3]; assign best_daq_lum to {"PB60:BDAQ:OWL","PB60:BDAQ:DAY","PB60:BDAQ:SWING"}; float last_tot_lum[3]; assign last_tot_lum to {"PB60:LTOT:OWL","PB60:LTOT:DAY","PB60:LTOT:SWING"}; monitor last_tot_lum; float last_col_lum[3]; assign last_col_lum to {"PB60:LCOL:OWL","PB60:LCOL:DAY","PB60:LCOL:SWING"}; monitor last_col_lum; float last_daq_lum[3]; assign last_daq_lum to {"PB60:LDAQ:OWL","PB60:LDAQ:DAY","PB60:LDAQ:SWING"}; monitor last_daq_lum; float last_3_shifts; assign last_3_shifts to "PB60:LCOL:LAST3"; monitor last_3_shifts; float shift_peak; assign shift_peak to "PB60:PEAK:SHIFT"; monitor shift_peak; float last_peaks[3]; assign last_peaks to {"PB60:PEAK:OWL","PB60:PEAK:DAY","PB60:PEAK:SWING"}; monitor last_peaks; short her_status; assign her_status to "HB60:STATE:ACT"; monitor her_status; short her_status_old; short ler_status; assign ler_status to "LB60:STATE:ACT"; monitor ler_status; short ler_status_old; short ler_logic; assign ler_logic to "LB60:INJECT:OVERRIDE"; short her_logic; assign her_logic to "HB60:INJECT:OVERRIDE"; short status_switch; float fill_thresh_l; assign fill_thresh_l to "LB60:BNCHCURR:FILL:THRESH"; monitor fill_thresh_l; float fill_thresh_h; assign fill_thresh_h to "HB60:BNCHCURR:FILL:THRESH"; monitor fill_thresh_h; float chisq_cut; assign chisq_cut to "PB60:LUM:LINEAR:REQ"; monitor chisq_cut; float chisq; assign chisq to "PB60:LUM:FITRVAL"; long manual_base; assign manual_base to "PB60:LUM:TIMEBASE:REQ"; monitor manual_base; long max_base; assign max_base to "PB60:LUM:TIMEBASE:MAX"; monitor max_base; long lum_lifemax; assign lum_lifemax to "PB60:LUM:LIFEMAX"; monitor lum_lifemax; long last_interval; assign last_interval to "PB60:LUM:TIMEBASE:ACT"; double lum_lifetime; assign lum_lifetime to "PB60:LUM:LIFETIME"; float current_slope; assign current_slope to "PB60:LUM:DLDT"; short beam_on_ler; assign beam_on_ler to "LB60:DCCT:BEAMON"; monitor beam_on_ler; short beam_on_her; assign beam_on_her to "HB60:DCCT:BEAMON"; monitor beam_on_her; short mpg_bit3_status; assign mpg_bit3_status to "PB60:MPG:BIT3STAT"; monitor mpg_bit3_status; short her_bit3_status; assign her_bit3_status to "HB60:BXBCM:BIT3STAT"; monitor her_bit3_status; short ler_bit3_status; assign ler_bit3_status to "LB60:BXBCM:BIT3STAT"; monitor ler_bit3_status; #define GOOD_BIT3_STATUS 16 /* "Always on" only */ /* Local variables */ long i; long j; long k; long lum_buf_index; long available_lum; %{ #define LUM_BUF_LEN 2000 typedef struct { double rel_time; double lum_value; }lumData_ts; /* File scope variables */ static lumData_ts lum_aps[LUM_BUF_LEN]; /* Array of structures */ TS_STAMP start_time; /* Time base */ TS_STAMP now; static long minute_of_day(TS_STAMP ptime) { struct tsDetail time_s; long i_min; /* long i_sec; */ tsStampToLocal(ptime, &time_s); i_min = time_s.hours*60+time_s.minutes; /* i_sec = time_s.seconds; printf("Bic_Lum MM:SS = %d:%d\n",i_min,i_sec); */ return(i_min); } /*============================================================================ Abs: Perform a linear fit to see if the decay is still flat Name: lum_linfit Type: int Args: CHISQ Use: Chi square of fit Type: Float Acc: Write Mech: By reference SLOPE Use: Slope of the fit line Type: float Acc: Write Mech: By reference INTERCEPT Use: Intercept of the fit line Type: float Acc: Write Mech: By reference CHISQ_CUT Use: Cutoff for CHISQ Type: Float Acc: Write Mech: By value NPTS Use: Number of points to fit Type: int Acc: read Mech: By value Rem: Butchered from REF_NAUTIL:BEV_LINFIT.FOR Side: Sets TIME_DIFF for GET_BASE_TIME Ret: Good or Bad fit (1 or 0) ============================================================================*/ int lum_linfit(float* chisq, float* slope, float* intercept, float chisq_cut, long int npts) { double sum=0, sumx = 0, sumx2 = 0, sumxy = 0, sumy = 0, sumy2 = 0; double xi=0, yi = 0; double delta,variance,a,b; double x_zero; /* Use x_zero in calculation below to get rid if problems after the display is up for a long time */ int i,j; *slope = *intercept = *chisq = 0.; if(beam_on_ler==0 || beam_on_her==0) return(0); if(npts<2) return(0); x_zero = lum_aps[lum_buf_index].rel_time; for(i=0, j=lum_buf_index; i=0.) *chisq = sqrt(variance); else return(0); if(*chisq<=chisq_cut) return(1); else return(0); } void get_base_time() { short j,fit_ok; float current_int; if (available_lum > 1) { /* First just age the last_interval and find the data associated with that age */ last_interval++; /* Check for a manual determination */ if(manual_base>0) last_interval=manual_base; /* Check the user-determined maximum time */ else if(max_base>0 && last_interval>max_base) last_interval=max_base; /* Check for enough data available */ if(last_interval>available_lum) last_interval=available_lum; fit_ok = lum_linfit(&chisq, ¤t_slope, ¤t_int, chisq_cut, last_interval); if(manual_base<=0 && fit_ok==0) for(j=last_interval/2; j>2&&last_interval>3; j/=2) { /* Adjust the fit interval up or down depending on the last fit results */ if (fit_ok==0) last_interval-=j; else last_interval+=j; fit_ok = lum_linfit(&chisq, ¤t_slope, ¤t_int, chisq_cut, last_interval); } } } }% ss lumSpinSS { state init { when(1) /* Start off by zeroing the luminosity array */ { %% taskDelay(300); for(i=0;i0.) num_saved_h++; pvPut(num_saved_h); epicsPrintf("HER NUMBUCK set to %i\n",num_saved_h); if(num_saved_h>10) { /* remove BLAST AWAY for more channels */ her_logic = 0; pvPut(her_logic); } } } state spin when(ler_status!=ler_status_old) { ler_status_old = ler_status; if(ler_status_old==Injecting) { status_switch = 1; pvGet(goal_l); num_saved_l=0; for(i=0;i0.) num_saved_l++; pvPut(num_saved_l); epicsPrintf("LER NUMBUCK set to %i\n",num_saved_l); if(num_saved_l>10) { /* remove BLAST AWAY for more channels */ ler_logic = 0; pvPut(ler_logic); } } } state spin when(status_switch==1) { status_switch=num_saved_o=0; for(i=0;i0. && goal_h[i]>0.) num_saved_o++; pvPut(num_saved_o); epicsPrintf("Goal overlap number set to %i \n",num_saved_o); } state spin when(ped_ler_old!=dcct_ler_ped) { if ( dcct_ler_ped > 10. || dcct_ler_ped < -10. ) { dcct_ler_ped = ped_ler_old; pvPut(dcct_ler_ped); PUBLISH_LER("Bad pedestal reset to old value") } ped_ler_old = dcct_ler_ped; } state spin when(ped_her_old!=dcct_her_ped) { if ( dcct_her_ped > 10. || dcct_her_ped < -10. ) { dcct_her_ped = ped_her_old; pvPut(dcct_her_ped); PUBLISH_HER("Bad pedestal reset to old value") } ped_her_old = dcct_her_ped; } state spin when(1) { /* Get the individual bucket current arrays */ pvGet(bc_l); pvGet(bc_h); now = pvTimeStamp(bc_h); /* Initialize the counters */ num_filled_l = num_filled_h = num_filled_o = 0; /* Loop over all the buckets, counting "full" buckets in each */ /* ring and "overlap" buckets, i.e. where both rings have */ /* current in the same bucket; also fill an array with I+*I- */ /* for identifying the luminosity distribution */ sq_i_l = 0.; sq_i_h = 0.; sq_i_hl = 0.; for(i=0;ifill_thresh_l); if(j) { num_filled_l++; sq_i_l += bc_l[i]*bc_l[i]; } if(bc_h[i]>fill_thresh_h) { num_filled_h++; sq_i_h += bc_h[i]*bc_h[i]; if(j) { num_filled_o++; bc_o[i]=bc_l[i]*bc_h[i]; sq_i_hl += bc_o[i]; } else bc_o[i]=0.; } else bc_o[i]=0.; }; pvPut(num_filled_l); pvPut(num_filled_h); pvPut(num_filled_o); pvPut(bc_o); /* Convert from microamperes squared to milliamperes squared */ sq_i_l *= .000001; pvPut(sq_i_l); sq_i_h *= .000001; pvPut(sq_i_h); sq_i_hl *= .000001; pvPut(sq_i_hl); /* Now calculate the "per-bucket" luminosity */ if(num_filled_o>0) buck_lum=dcct_l*dcct_h/num_filled_o; else buck_lum = 0.; pvPut(buck_lum); /* Get the current minute */ %% i_min = minute_of_day(now); splumval_inv = dcct_l*dcct_h*num_saved_o; if(splumval_inv > 10. && lumval > 0.) splumval=lumval*num_saved_h*num_saved_l/splumval_inv; else splumval=0.; if(splumval>5.) splumval=5.; pvPut(splumval); if(lumval>shift_peak) { shift_peak=lumval; pvPut(shift_peak); } if(i_min!=i_min_old) { i_min_old = i_min; /* printf("%ld minutes since midnight\n",i_min); */ day_lum[i_min]=lumval; day_heri[i_min]=dcct_h; day_leri[i_min]=dcct_l; /* day_splum[i_min]=splumval; */ /* pvPut(day_splum); */ pvPut(day_leri); pvPut(day_heri); pvPut(day_lum); %% memmove((void*)HERI_LOC,day_heri, sizeof(float)*N_MINUTES); %% memmove((void*)LERI_LOC,day_leri, sizeof(float)*N_MINUTES); %% memmove((void*)LUM_LOC, day_lum, sizeof(float)*N_MINUTES); for(i=0;ibest_tot_lum[old_shift]) { best_tot_lum[old_shift]=current_tot_lum; pvPut(best_tot_lum[old_shift]); } last_tot_lum[old_shift]=current_tot_lum; pvPut(last_tot_lum[old_shift]); pvGet(best_col_lum[old_shift]); if(current_col_lum>best_col_lum[old_shift]) { best_col_lum[old_shift]=current_col_lum; pvPut(best_col_lum[old_shift]); } last_col_lum[old_shift]=current_col_lum; pvPut(last_col_lum[old_shift]); pvGet(best_daq_lum[old_shift]); if(current_daq_lum>best_daq_lum[old_shift]) { best_daq_lum[old_shift]=current_daq_lum; pvPut(best_daq_lum[old_shift]); } last_daq_lum[old_shift]=current_daq_lum; pvPut(last_daq_lum[old_shift]); last_peaks[old_shift] = shift_peak; pvPut(last_peaks[old_shift]); old_shift=current_shift; current_tot_lum=0; current_col_lum=0; pvPut(current_col_lum); current_daq_lum=0; pvPut(current_daq_lum); shift_peak = lumval; pvPut(shift_peak); last_3_shifts=0; for(i=0;i<3;i++) last_3_shifts+=last_col_lum[i]; pvPut(last_3_shifts); } pvPut(current_shift); /* Normal minute update */ current_tot_lum+=lumval*INT_CONST; pvPut(current_tot_lum); current_daq_lum+=lumdaq*INT_CONST; pvPut(current_daq_lum); current_col_lum+=lumcoll*INT_CONST; pvPut(current_col_lum); if(mpg_bit3_status!=GOOD_BIT3_STATUS) { PUBLISH_LER("MPG Bit-3 memory problem"); PUBLISH_HER("MPG Bit-3 memory problem"); } if(her_bit3_status!=GOOD_BIT3_STATUS) { PUBLISH_HER("HER Bit-3 memory problem"); } if(ler_bit3_status!=GOOD_BIT3_STATUS) { PUBLISH_LER("LER Bit-3 memory problem"); } } /* end of new minute reached */ /* And finally, wait 1 second before the next loop */ lum_buf_index++; if(available_lum=LUM_BUF_LEN) lum_buf_index=0; TsDiffAsDouble(&lum_aps[lum_buf_index].rel_time, &now,&start_time); if(lum_aps[lum_buf_index].rel_time>3.E7) { start_time=now; lum_aps[lum_buf_index].rel_time=0.; } /* lum_aps[lum_buf_index].lum_value = lumval; */ lum_aps[lum_buf_index].lum_value = lummed; get_base_time(); pvPut(current_slope); if(fabs(current_slope) > .001) lum_lifetime= -.01667*lum_aps[lum_buf_index].lum_value/current_slope; else lum_lifetime = 0.; if(lum_lifetime > lum_lifemax) lum_lifetime = lum_lifemax; if(lum_lifetime < 0.) lum_lifetime=0.; pvPut(last_interval); pvPut(chisq); pvPut(lum_lifetime); %% taskDelay(60); /* Pick up i_min in case it was set to 2000 or 3000 for a reset */ pvGet(i_min); } state spin } } exit{}