Compare revisions

62f0f2e4 · 62f0f2e4 · 62f0f2e4 · 62f0f2e4 · 62f0f2e4 · 62f0f2e4
--- a/docs/tutorials/data/DLCs_onshore/Main.txt
+++ b/docs/tutorials/data/DLCs_onshore/Main.txt
+																																					
+Master spreadsheets to generate the set of spreadsheets required as inputs to the DLB calculator.																																					
+"Each sheet defines the tags of a DLC, except the main one. The main sheet defines: wind turbine parameters, default tags values, and gusts and turbulence definitions."																																					
+"Tags are devided into 3 categories: constants (C), variables (V), and functions (F).  The category is specified in the line above the tag."																																					
+Constants do not change in a DLC. Variables define the number of cases within a DLC through their combinations. Functions are tags that depends on other tags through and expression. 																																					
+																																					
+Parameters:	Vrate	Vout																																			
+	12	26																																			
+																																					
+Default constants:	[ref_ti]	[ref_wind_speed]	[tsr]	[hub_height]	[diameter]	[t0]	[wdir]	[shear_exp]	[out_format]	[gust]	[gust_type]	[G_A]	[G_phi0]	[G_t0]	[G_T]	[Rotor azimuth]	[Free shaft rot]	[init_wr]	[Pitch 1 DLC22b]	[Rotor locked]	[Time stuck DLC22b]	[Cut-in time]	[Stop type]	[Pitvel 1]	[Pitvel 2]	[Grid loss time]	[Time pitch runaway]	[Induction]	[Dyn stall]	[dis_setbeta]	[long_scale_param]	[t_flap_on]	[turb_format]	[staircase]	[Rotor azimuth]	[sim_time]	[Cut-out time]
+	0.16	50	8.0	90	178	100	0	0.2	hawc_binary	;						0		0.5	0	;	-1	1	1	4	6	10000	10000	1	2	 	42	20	1	;	0	600	10000
+																																					
+																																					
+Default functions:	[Turb base name]	[time stop]	[turb_dx]	[wsp factor]	[wind_ramp_t1]	[wind_ramp_factor1]	[time_start]																														
+	"""turb_wsp[wsp]_s[seed]"""	[t0]+[sim_time]	"[wsp]*[sim_time]/8192,0"	[tsr]/[wsp]	[t0]	[wsp factor]	[t0]																														
+																																					
+																																					
+																																					
+Gusts:																																					
+EOG	"min([1,35*(0,8*1,4*[ref_wind_speed]-[wsp]);3,3*[TI]*[wsp]/(1+0,1*[diameter]/[long_scale_param])])"																																				
+ECD	15																																				
+EWS	"(2,5+0,2*6,4*[TI]*[wsp]*([diameter]/[long_scale_param])**0,25)/[diameter]"																																				
+EDC	"4*arctan([TI]/(1+0,1*[diameter]/[long_scale_param]))*180/pi"																																				
+																																					
+Turbulence:																																					
+NTM	"([ref_ti]*(0,75*[wsp]+5,6))/[wsp]"																																				
+ETM	"2*[ref_ti]*(0,072*(0,2*[ref_wind_speed]/2+3)*([wsp]/2-4)+10)/[wsp]"																																				
+																																					
+	Wind speeds:																																				
+	4																																				
+	6																																				
+	8																																				
+	10																																				
+	12																																				
+	14																																				
+	16																																				
+	18																																				
+	20																																				
+	22																																				
+	24																																				
+	26																																				
+	28																																				
+	30																																				
+	32																																				
+	34																																				
+	36																																				
+	38																																				
+	40																																				
+	42																																				
+	44																																				
+	46																																				
+	48																																				
+	50																																				
--- a/docs/tutorials/data/dtu10mw_master_onshore.htc
+++ b/docs/tutorials/data/dtu10mw_master_onshore.htc
+; DTU_10MW_RWT, version 5, 04-21-15, anyd 
+;
+; Demo master file for Wind Energy Toolbox
+; 
+begin simulation;
+  time_stop    [time stop];
+  solvertype   1 ;    (newmark)
+  on_no_convergence continue ;
+  convergence_limits 1E3 1.0 1E-7 ;
+  logfile ./logfiles/[Case folder]/[Case id.].log ;
+  begin newmark;
+    deltat    0.02;  
+  end newmark;
+end simulation;
+;
+;---------------------------------------------------------------------------------------------------------------------------------------------------------------- 
+begin new_htc_structure;
+;--------------------------------------------------------------------------------------------------
+[staircase]  beam_output_file_name  ./[eigenfreq_dir][Case folder]/[Case id.]/[Case id.]_beam.dat;
+[staircase]  body_output_file_name  ./[eigenfreq_dir][Case folder]/[Case id.]/[Case id.]_body.dat;
+[staircase]  struct_inertia_output_file_name ./[eigenfreq_dir][Case folder]/[Case id.]/[Case id.]_struct.dat;
+[staircase]  body_eigenanalysis_file_name ./[eigenfreq_dir][Case folder]/[Case id.]/[Case id.]_body_eigen.dat;
+[staircase]  structure_eigenanalysis_file_name ./[eigenfreq_dir][Case folder]/[Case id.]/[Case id.]_strc_eigen.dat;
+;---------------------------------------------------------------------------------------------------
+  begin main_body;         tower 115m
+    name        tower ;            
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution     c2_def ;
+    damping_posdef   0.0 0.0 0.0 4.12E-03 4.12E-03 4.5E-04  ; Mx My Mz Kx Ky Kz , Ms raises overall level, Ks raises high freguency level "tuned by Larh"
+     begin timoschenko_input;
+      filename ./data/DTU_10MW_RWT_Tower_st.dat;
+      set 1 1 ; 
+    end timoschenko_input;
+    begin c2_def;              Definition of centerline (main_body coordinates)
+      nsec 11;
+      sec	1	0	0	0.00	0	;  x,y,z,twist
+      sec	2	0	0	-11.50	0	;
+      sec	3	0	0	-23.00	0	;
+      sec	4	0	0	-34.50	0	;
+      sec	5	0	0	-46.00	0	;
+      sec	6	0	0	-57.50	0	;
+      sec	7	0	0	-69.00	0	;
+      sec	8	0	0	-80.50	0	;
+      sec	9	0	0	-92.00	0	 ;
+      sec	10	0	0	-103.50	0	;
+      sec	11	0	0	-115.63	0	;	  
+     end c2_def ;
+    end main_body;
+;
+  begin main_body;
+    name        towertop ;              
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution     c2_def ;
+    damping_posdef  0.0  0.0  0.0  7.00E-03  7.00E-03  7.00E-03  ;   "changed by Larh"	
+	concentrated_mass	2.0	0.0	2.6870E+00	3.0061E-01	4.4604E+05	4.1060E+06	4.1060E+05	4.1060E+06	;	Nacelle mass and inertia "corrected by Anyd 25/4/13"
+	begin timoschenko_input;
+      filename ./data/DTU_10MW_RWT_Towertop_st.dat ;
+      set 1 2 ;                
+    end timoschenko_input;
+    begin c2_def;              Definition of centerline (main_body coordinates)
+      nsec 2;
+      sec 1 0.0 0.0  0.0    0.0 ; x,y,z,twist
+      sec 2 0.0 0.0 -2.75   0.0 ; 
+    end c2_def ;
+  end main_body;
+;
+  begin main_body;
+    name        shaft ;              
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution     c2_def ;
+	damping_posdef  0.0 0.0 0.0 4.65E-04  4.65E-04  3.983E-03 ; "tuned by Anyd 23/5/13 to 31.45 log decr. damping for free free with stiff rotor and tower"
+    concentrated_mass	1.0	0.0	0.0	0.0	0.0	0.0	0.0	3.751E+06	;	generator equivalent slow shaft "re_tuned by Anyd 20/2/13"  
+    concentrated_mass	5.0	0.0	0.0	0.0	1.0552E+05	0.0	0.0	3.257E+05	;	hub mass and inertia;	"re_tuned by Anyd 20/2/13"  
+	begin timoschenko_input;
+      filename ./data/DTU_10MW_RWT_Shaft_st.dat ;
+      set 1 1 ;                
+    end timoschenko_input;
+    begin c2_def;              Definition of centerline (main_body coordinates)
+      nsec 5;
+      sec 1 0.0 0.0 0.0 0.0 ; Tower top x,y,z,twist
+      sec 2 0.0 0.0 1.5 0.0 ; 
+      sec 3 0.0 0.0 3.0 0.0 ; 
+      sec 4 0.0 0.0 4.4	0.0 ; Main bearing
+      sec 5 0.0 0.0 7.1 0.0 ; Rotor centre
+    end c2_def ;
+  end main_body;	
+;
+  begin main_body;
+    name        hub1 ;              
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution     c2_def ;
+    damping_posdef  0.0  0.0  0.0  3.00E-06  3.00E-06  2.00E-05;   "changed by Larh"	
+	begin timoschenko_input;
+      filename ./data/DTU_10MW_RWT_Hub_st.dat ;
+      set 1 2 ;                
+    end timoschenko_input;
+    begin c2_def;              Definition of centerline (main_body coordinates)
+      nsec 2;
+      sec 1 0.0 0.0 0.0 0.0 ; x,y,z,twist
+      sec 2 0.0 0.0 2.8 0.0 ; 
+    end c2_def ;
+  end main_body;
+;
+  begin main_body;
+    name           hub2 ;
+    copy_main_body hub1;
+  end main_body;
+;
+  begin main_body;
+    name           hub3 ;
+    copy_main_body hub1 ;
+  end main_body;
+;
+  begin main_body;
+    name        blade1 ;        
+    type        timoschenko ;
+    nbodies     10 ;
+    node_distribution    c2_def;
+	damping_posdef   0.0 0.0 0.0 1.53e-3 2.55e-3 3.3e-4 ; " 3% damping tuned by tkim 23/03/13 unable to fit 3rd and higher mode"
+    begin timoschenko_input ;
+      filename ./data/DTU_10MW_RWT_Blade_st.dat;
+      set 1 1 ;                set subset
+    end timoschenko_input;
+    begin c2_def;              Definition of centerline (main_body coordinates)
+		nsec 27 ;
+			sec	1	0.00000E+00	7.00600E-05	4.44089E-16	-1.45000E+01	;
+			sec	2	-2.06477E-05	-1.22119E-02	3.00000E+00	-1.45000E+01	;
+			sec	3	-7.28810E-03	-2.49251E-02	6.00000E+00	-1.44851E+01	;
+			sec	4	-1.89235E-02	-2.73351E-02	7.00004E+00	-1.44610E+01	;
+			sec	5	-5.41282E-02	-2.82163E-02	8.70051E+00	-1.43388E+01	;
+			sec	6	-1.26633E-01	-2.13210E-02	1.04020E+01	-1.40201E+01	;
+			sec	7	-2.25666E-01	-1.28378E-02	1.22046E+01	-1.33904E+01	;
+			sec	8	-2.88563E-01	-7.70659E-03	1.32065E+01	-1.29371E+01	;
+			sec	9	-3.99194E-01	-4.88317E-03	1.50100E+01	-1.19445E+01	;
+			sec	10	-5.76634E-01	-1.80296E-02	1.82151E+01	-9.98243E+00	;
+			sec	11	-7.07136E-01	-5.01772E-02	2.14178E+01	-8.45147E+00	;
+			sec	12	-7.91081E-01	-9.41228E-02	2.46189E+01	-7.46417E+00	;
+			sec	13	-8.37195E-01	-1.48880E-01	2.78193E+01	-6.72916E+00	;
+			sec	14	-8.53948E-01	-2.14514E-01	3.10194E+01	-6.08842E+00	;
+			sec	15	-8.49367E-01	-2.90618E-01	3.42197E+01	-5.49322E+00	;
+			sec	16	-7.93920E-01	-4.62574E-01	4.02204E+01	-4.39222E+00	;
+			sec	17	-7.16284E-01	-6.88437E-01	4.66217E+01	-3.09315E+00	;
+			sec	18	-6.34358E-01	-9.60017E-01	5.30232E+01	-1.75629E+00	;
+			sec	19	-5.53179E-01	-1.28424E+00	5.94245E+01	-5.00650E-01	;
+			sec	20	-4.75422E-01	-1.66402E+00	6.58255E+01	6.01964E-01	;
+			sec	21	-4.03180E-01	-2.10743E+00	7.22261E+01	1.55560E+00	;
+			sec	22	-3.30085E-01	-2.65630E+00	7.90266E+01	2.51935E+00	;
+			sec	23	-3.10140E-01	-2.78882E+00	8.05267E+01	2.72950E+00	;
+			sec	24	-2.86719E-01	-2.92517E+00	8.20271E+01	2.93201E+00	;
+			sec	25	-2.55823E-01	-3.06577E+00	8.35274E+01	3.11874E+00	;
+			sec	26	-2.07891E-01	-3.20952E+00	8.50277E+01	3.28847E+00	;
+			sec	27	-8.98940E-02	-3.33685E+00	8.63655E+01	3.42796E+00	;
+     end c2_def ;                                 
+   end main_body;
+;
+  begin main_body;
+    name           blade2 ;
+    copy_main_body blade1;
+  end main_body;
+;
+  begin main_body;
+    name           blade3 ;
+    copy_main_body blade1 ;
+  end main_body;
+;-------------------------------------------------------------------------------------------------------------------------------                   
+;
+  begin orientation;
+    begin base;
+      body   tower;
+      inipos        0.0 0.0 0.0 ;         initial position of node 1
+      body_eulerang 0.0 0.0 0.0;
+    end base;
+; 
+    begin relative;
+      body1  tower last;
+      body2  towertop 1;
+      body2_eulerang 0.0 0.0 0.0; 
+    end relative;
+;
+    begin relative;
+      body1  towertop last;
+      body2  shaft 1;
+      body2_eulerang 90.0 0.0 0.0; 
+      body2_eulerang 5.0 0.0 0.0;    5 deg tilt angle
+      body2_eulerang 0.0 0.0 [Rotor azimuth];
+[Free shaft rot]      mbdy2_ini_rotvec_d1 0.0 0.0 -1.0 [init_wr] ;
+    end relative;
+;
+    begin relative;
+      body1  shaft last;         
+      body2  hub1 1;
+      body2_eulerang -90.0 0.0 0.0;    
+      body2_eulerang 0.0 180.0 0.0;    
+      body2_eulerang 2.5 0.0 0.0;      2.5deg cone angle
+    end relative;
+;
+    begin relative;
+      body1  shaft last;         
+      body2  hub2 1;
+      body2_eulerang -90.0 0.0 0.0;    
+      body2_eulerang 0.0 60.0 0.0;   
+      body2_eulerang 2.5 0.0 0.0;      2.5deg cone angle
+    end relative;
+;
+    begin relative;
+      body1  shaft last;         
+      body2  hub3 1;
+      body2_eulerang -90.0 0.0 0.0;    
+      body2_eulerang 0.0 -60.0 0.0;    
+      body2_eulerang 2.5 0.0 0.0;      2.5deg cone angle
+    end relative;
+;
+    begin relative;
+      body1  hub1 last;         
+      body2  blade1 1;
+      body2_eulerang 0.0 0.0 0;    
+    end relative;
+;
+    begin relative;
+      body1  hub2 last;         
+      body2  blade2 1;
+      body2_eulerang 0.0 0.0 0.0;    
+    end relative;
+;
+    begin relative;
+      body1  hub3 last;         
+      body2  blade3 1;
+      body2_eulerang 0.0 0.0 0.0;    
+    end relative;
+;
+ 	end orientation;
+;-------------------------------------------------------------------------------------------------------------------------------
+begin constraint;   
+;
+    begin fix0;  fixed to ground in translation and rotation of node 1
+      body tower;
+    end fix0;
+;
+     begin fix1;
+		   body1 tower last ;
+		   body2 towertop 1;
+		 end fix1;
+;
+[Free shaft rot]    begin bearing1;                       free bearing
+[Free shaft rot]      name  shaft_rot;
+[Free shaft rot]      body1 towertop last;
+[Free shaft rot]      body2 shaft 1;
+[Free shaft rot]      bearing_vector 2 0.0 0.0 -1.0;        x=coo (0=global.1=body1.2=body2) vector in body2 coordinates where the free rotation is present
+[Free shaft rot]   end bearing1; 
+;
+[Rotor locked]    begin bearing3;                       free bearing
+[Rotor locked]      name  shaft_rot;
+[Rotor locked]      body1 towertop last;
+[Rotor locked]      body2 shaft 1;
+[Rotor locked]      bearing_vector 2 0.0 0.0 -1.0;        x=coo (0=global.1=body1.2=body2) vector in body2 coordinates where the free rotation is present
+[Rotor locked]      omegas 0.0 ;
+[Rotor locked]   end bearing3; 
+;
+     begin fix1;
+		   body1 shaft last ;
+		   body2 hub1 1;
+		 end fix1;
+;
+     begin fix1;
+		   body1 shaft last ;
+		   body2 hub2 1;
+		 end fix1;
+;
+     begin fix1;
+		   body1 shaft last ;
+		   body2 hub3 1;
+		 end fix1; 
+;	
+    begin bearing2;
+      name pitch1;		
+      body1 hub1 last;
+     body2 blade1 1;
+			bearing_vector 2 0.0 0.0 -1.0;
+   end bearing2;
+;
+    begin bearing2;
+      name pitch2;		
+      body1 hub2 last;
+      body2 blade2 1;
+			bearing_vector 2 0.0 0.0 -1.0;
+    end bearing2;
+;
+    begin bearing2;
+      name pitch3;		
+      body1 hub3 last;
+      body2 blade3 1;
+			bearing_vector 2 0.0 0.0 -1.0;
+    end bearing2;
+end constraint;
+;
+end new_htc_structure;
+;---------------------------------------------------------------------------------------------------------------------------------------------------------------- 
+begin wind ;
+  density                 1.225 ;
+  wsp                     [wsp]   ;
+  tint                    [TI]   ;
+  horizontal_input        1     ;
+  windfield_rotations     [wdir] 0.0 0.0 ;    yaw, tilt, rotation
+  center_pos0             0.0 0.0 -119 ; hub heigth
+  shear_format            3 [shear_exp] ;
+  turb_format             [turb_format]     ;  0=none, 1=mann,2=flex
+  tower_shadow_method     3     ;  0=none, 1=potential flow, 2=jet
+  scale_time_start       [t0] ;						   
+  wind_ramp_factor        0.0 [t0] [wsp factor] 1.0 ;
+ [gust] iec_gust               [gust_type] [G_A] [G_phi0] [G_t0] [G_T] ;
+; 
+[staircase]  wind_ramp_abs  400.0  401.0  0.0   1.0 ;   wsp. after the step:  5.0 
+[staircase]  wind_ramp_abs  501.0  502.0  0.0   1.0 ;   wsp. after the step:  6.0 
+[staircase]  wind_ramp_abs  602.0  603.0  0.0   1.0 ;   wsp. after the step:  7.0 
+[staircase]  wind_ramp_abs  703.0  704.0  0.0   1.0 ;   wsp. after the step:  8.0 
+[staircase]  wind_ramp_abs  804.0  805.0  0.0   1.0 ;   wsp. after the step:  9.0 
+[staircase]  wind_ramp_abs  905.0  906.0  0.0   1.0 ;   wsp. after the step: 10.0 
+[staircase]  wind_ramp_abs  1006.0  1007.0  0.0   1.0 ;   wsp. after the step: 11.0 
+[staircase]  wind_ramp_abs  1107.0  1108.0  0.0   1.0 ;   wsp. after the step: 12.0 
+[staircase]  wind_ramp_abs  1208.0  1209.0  0.0   1.0 ;   wsp. after the step: 13.0 
+[staircase]  wind_ramp_abs  1309.0  1310.0  0.0   1.0 ;   wsp. after the step: 14.0 
+[staircase]  wind_ramp_abs  1410.0  1411.0  0.0   1.0 ;   wsp. after the step: 15.0 
+[staircase]  wind_ramp_abs  1511.0  1512.0  0.0   1.0 ;   wsp. after the step: 16.0 
+[staircase]  wind_ramp_abs  1612.0  1613.0  0.0   1.0 ;   wsp. after the step: 17.0 
+[staircase]  wind_ramp_abs  1713.0  1714.0  0.0   1.0 ;   wsp. after the step: 18.0 
+[staircase]  wind_ramp_abs  1814.0  1815.0  0.0   1.0 ;   wsp. after the step: 19.0 
+[staircase]  wind_ramp_abs  1915.0  1916.0  0.0   1.0 ;   wsp. after the step: 20.0 
+[staircase]  wind_ramp_abs  2016.0  2017.0  0.0   1.0 ;   wsp. after the step: 21.0 
+[staircase]  wind_ramp_abs  2117.0  2118.0  0.0   1.0 ;   wsp. after the step: 22.0 
+[staircase]  wind_ramp_abs  2218.0  2219.0  0.0   1.0 ;   wsp. after the step: 23.0 
+[staircase]  wind_ramp_abs  2319.0  2320.0  0.0   1.0 ;   wsp. after the step: 24.0 
+[staircase]  wind_ramp_abs  2420.0  2421.0  0.0   1.0 ;   wsp. after the step: 25.0 
+  ;
+  begin mann ;
+    create_turb_parameters 29.4 1.0 3.9 [seed] 1.0 ;      L, alfaeps, gamma, seed, highfrq compensation
+    filename_u    ./[turb_dir][Turb base name]u.bin ;
+    filename_v    ./[turb_dir][Turb base name]v.bin ;
+    filename_w    ./[turb_dir][Turb base name]w.bin ;
+    box_dim_u    8192 [turb_dx] ;
+    box_dim_v    32 6.5;
+    box_dim_w    32 6.5;
+    std_scaling   1.0 0.7 0.5 ;
+  end mann ;
+;  
+  begin tower_shadow_potential_2;
+    tower_mbdy_link tower;
+    nsec  2;
+    radius      0.0  4.15 ;
+    radius     115.63 2.75 ;
+  end tower_shadow_potential_2;
+end wind;
+;
+begin aerodrag ;
+  begin aerodrag_element ;
+    mbdy_name tower;
+    aerodrag_sections uniform 10 ;
+    nsec 2 ;
+    sec 0.0 0.6 8.3 ;  tower bottom
+    sec 115.63 0.6 5.5 ;  tower top
+  end aerodrag_element;
+;
+  begin aerodrag_element ;        Nacelle drag side
+    mbdy_name shaft;
+    aerodrag_sections uniform 2 ;
+    nsec 2 ;
+    sec 0.0   0.8 10.0 ;  
+    sec 7.01  0.8 10.0 ;  
+  end aerodrag_element;
+end aerodrag;
+;
+begin aero ;
+  nblades  3;
+  hub_vec shaft -3 ;         rotor rotation vector (normally shaft composant directed from pressure to sustion side)
+  link 1 mbdy_c2_def blade1;
+  link 2 mbdy_c2_def blade2;
+  link 3 mbdy_c2_def blade3;
+  ae_filename        ./data/DTU_10MW_RWT_ae.dat ;
+  pc_filename        ./data/DTU_10MW_RWT_pc.dat ;
+  induction_method   [Induction] ;     0=none, 1=normal
+  aerocalc_method    1 ;     0=ingen aerodynamic, 1=med aerodynamic
+  aerosections       50 ; def. 50
+  ae_sets            1 1 1;
+  tiploss_method     1 ;     0=none, 1=prandtl
+  dynstall_method    [Dyn stall];     0=none, 1=stig øye method,2=mhh method
+;  
+end aero ;
+;-------------------------------------------------------------------------------------------------  
+begin dll;
+;
+  begin type2_dll;
+    name dtu_we_controller ;
+    filename  ./control/dtu_we_controller.dll ;
+    dll_subroutine_init init_regulation ;
+    dll_subroutine_update update_regulation ;
+    arraysizes_init  100 1 ;
+    arraysizes_update  50 100 ;
+	begin init ;
+       ; Overall parameters
+      constant   1 10000.0    	; Rated power [kW]                         
+      constant   2   0.628    	; Minimum rotor speed [rad/s]
+      constant   3   1.005    	; Rated rotor speed [rad/s]
+      constant   4  15.6E+06  	; Maximum allowable generator torque [Nm]
+      constant   5  100.0    	; Minimum pitch angle, theta_min [deg], 
+								; if |theta_min|>90, then a table of <wsp,theta_min> is read ;
+								; from a file named 'wptable.n', where n=int(theta_min)
+      constant   6  82.0    	; Maximum pitch angle [deg]
+      constant   7  10.0    	; Maximum pitch velocity operation [deg/s]
+      constant   8   0.4    	; Frequency of generator speed filter [Hz]
+      constant   9   0.7    	; Damping ratio of speed filter [-]
+      constant  10   1.92   	; Frequency of free-free DT torsion mode [Hz], if zero no notch filter used
+      ; Partial load control parameters
+      constant  11   1.175E+07  ; Optimal Cp tracking K factor [Nm/(rad/s)^2], ;
+                                ; Qg=K*Omega^2, K=eta*0.5*rho*A*Cp_opt*R^3/lambda_opt^3                     
+      constant  12   7.084E+07  ; Proportional gain of torque controller [Nm/(rad/s)]
+      constant  13   1.590E+07  ; Integral gain of torque controller [Nm/rad]
+      constant  14   0.0    	; Differential gain of torque controller [Nm/(rad/s^2)]
+;     Full load control parameters
+      constant  15   2      	; Generator control switch [1=constant power, 2=constant torque]
+      constant  16   1.304E+00  ; Proportional gain of pitch controller [rad/(rad/s)]
+      constant  17   3.511E-01  ; Integral gain of pitch controller [rad/rad]
+      constant  18   0.0    	; Differential gain of pitch controller [rad/(rad/s^2)]
+      constant  19   0.4e-8 	; Proportional power error gain [rad/W]
+      constant  20   0.4e-8 	; Integral power error gain [rad/(Ws)]
+      constant  21 	 1.135E+01  ; Coefficient of linear term in aerodynamic gain scheduling, KK1 [deg]
+      constant  22 	 4.007E+02  ; Coefficient of quadratic term in aerodynamic gain scheduling, KK2 [deg^2] &
+								; (if zero, KK1 = pitch angle at double gain)
+      constant  23   1.3    	; Relative speed for double nonlinear gain [-]
+;     Cut-in simulation parameters
+      constant  24  [Cut-in time] 		; Cut-in time [s]
+      constant  25  1.0   		; Time delay for soft start of torque [1/1P]
+;     Cut-out simulation parameters
+      constant  26  [Cut-out time]   ; Cut-out time [s]
+      constant  27   5.0 	; Time constant for linear torque cut-out [s]
+      constant  28  [Stop type]   	; Stop type [1=normal, 2=emergency]
+      constant  29  1.0  	; Time delay for pitch stop after shut-down signal [s]
+      constant  30  [Pitvel 1]  	; Maximum pitch velocity during initial period of stop [deg/s]
+      constant  31  3.0  	; Time period of initial pitch stop phase [s] (maintains pitch speed specified in constant 30)
+      constant  32  [Pitvel 2] 	; Maximum pitch velocity during final phase of stop [deg/s]
+;     Expert parameters (keep default values unless otherwise given)
+      constant  33   2.0  	; Lower angle above lowest minimum pitch angle for switch [deg]
+      constant  34   2.0  	; Upper angle above lowest minimum pitch angle for switch [deg], if equal then hard switch
+      constant  35  95.0  	; Ratio between filtered speed and reference speed for fully open torque limits [%]
+      constant  36   2.0  	; Time constant of 1st order filter on wind speed used for minimum pitch [1/1P]
+      constant  37   1.0  	; Time constant of 1st order filter on pitch angle used for gain scheduling [1/1P]
+;     Drivetrain damper
+      constant  38   0.0  	; Proportional gain of active DT damper [Nm/(rad/s)], requires frequency in input 10
+;	  Over speed
+	  constant  39  25.0  	; Overspeed percentage before initiating turbine controller alarm (shut-down) [%]
+;     Additional non-linear pitch control term (not used when all zero)
+	  constant  40   0.0  	; Err0 [rad/s] 
+	  constant  41   0.0  	; ErrDot0 [rad/s^2]
+	  constant  42   0.0  	; PitNonLin1 [rad/s]
+;     Storm control command
+	  constant 43   28.0  	; Wind speed 'Vstorm' above which derating of rotor speed is used [m/s]
+	  constant 44   28.0  	; Cut-out wind speed (only used for derating of rotor speed in storm) [m/s]	  
+;     Safety system parameters
+	  constant 45   30.0  ; Overspeed percentage before initiating safety system alarm (shut-down) [%]
+	  constant 46    1.5  ; Max low-pass filtered tower top acceleration level [m/s^2]
+;     Turbine parameter
+	  constant 47  178.0  ; Nominal rotor diameter [m]
+;     Parameters for rotor inertia reduction in variable speed region
+      constant 48    0.0  ; Proportional gain on rotor acceleration in variable speed region [Nm/(rad/s^2)] (not used when zero)
+;     Parameters for alternative partial load controller with PI regulated TSR tracking
+      constant 49    0.0  ; Optimal tip speed ratio [-] (only used when K=constant 11 = 0 otherwise  Qg=K*Omega^2 is used)
+;     Parameters for adding aerodynamic drivetrain damping on gain scheduling
+      constant 50    0.0  ; Proportional gain of aerodynamic DT damping [Nm/(rad/s)]
+      constant 51    0.0  ; Coefficient of linear term in aerodynamic DT damping scheduling, KK1 [deg]
+      constant 52    0.0  ; Coefficient of quadratic term in aerodynamic DT damping scheduling, KK2 [deg^2]
+	end init ;
+;
+    begin output ;
+      general time ; [s]     
+      constraint bearing1 shaft_rot 1 only 2 ; Drivetrain speed [rad/s]
+      constraint bearing2 pitch1 1 only 1; [rad]         
+      constraint bearing2 pitch2 1 only 1; [rad]                               
+      constraint bearing2 pitch3 1 only 1; [rad]                               
+      wind free_wind 1 0.0 0.0 -119 	 ; Global coordinates at hub height
+	  dll inpvec 2 2 					 ; Elec. power from generator servo .dll
+	  dll inpvec 2 8 					 ; Grid state flag from generator servo .dll
+	  mbdy state acc tower 10 1.0 global only 1 ; Tower top x-acceleration [m/s^2]
+	  mbdy state acc tower 10 1.0 global only 2 ; Tower top y-acceleration [m/s^2]
+    end output;     
+  end type2_dll;
+;
+   begin type2_dll;
+     name generator_servo ;
+     filename  ./control/generator_servo.dll ;
+     dll_subroutine_init init_generator_servo ;
+     dll_subroutine_update update_generator_servo ;
+     arraysizes_init  7 1 ;
+     arraysizes_update  4 8 ;
+ 	begin init ;
+       constant 1  20.0    ; Frequency of 2nd order servo model of generator-converter system [Hz]   
+       constant 2  0.9     ; Damping ratio 2nd order servo model of generator-converter system [-]
+       constant 3 15.6E+06 ; Maximum allowable LSS torque (pull-out torque) [Nm]
+       constant 4 0.94     ; Generator efficiency [-]
+       constant 5 1.0      ; Gearratio [-]
+       constant 6 0.0      ; Time for half value in softstart of torque [s]
+       constant 7 [Grid loss time]     ; Time for grid loss [s]
+     end init ;
+;
+     begin output;
+       general time                          	; Time [s]    
+       dll inpvec 1 1                        	; Electrical torque reference [Nm]  
+       constraint bearing1 shaft_rot 1 only 2	; Generator LSS speed [rad/s]   
+       mbdy momentvec shaft 1 1 shaft only 3 	; Shaft moment [kNm] (Qshaft)
+     end output;
+;
+     begin actions;    
+        mbdy moment_int shaft 1 -3 shaft towertop 2 ;   Generator LSS torque [Nm]
+     end actions;
+   end type2_dll;
+;
+   begin type2_dll;
+     name mech_brake ;
+     filename  ./control/mech_brake.dll ;
+     dll_subroutine_init init_mech_brake ;
+     dll_subroutine_update update_mech_brake ;
+     arraysizes_init    7 1 ;
+     arraysizes_update  4 6 ;
+ 	begin init ;
+      constant 1 2727252.0   ; Fully deployed maximum brake torque [Nm]
+      constant 2     100.0   ; Parameter alpha used in Q = tanh(omega*alpha), typically 1e2/Omega_nom
+      constant 3       0.625 ; Delay time for before brake starts to deploy [s] - from 5MW*1P_5/1P_10
+      constant 4       0.75  ; Time for brake to become fully deployed [s]
+     end init ;
+;
+     begin output;
+	   general time                           ; Time [s]
+	   constraint bearing1 shaft_rot 1 only 2 ; Generator LSS speed [rad/s]
+	   dll inpvec 1 25                        ; Command to deploy mechanical disc brake [0,1]
+     end output;
+;
+     begin actions;    
+        mbdy moment_int shaft 1 3 shaft towertop 2 ;   Brake LSS torque [Nm]
+     end actions;
+   end type2_dll;
+;
+  begin type2_dll;
+    name servo_with_limits ;
+    filename  ./control/servo_with_limits.dll ;
+    dll_subroutine_init init_servo_with_limits ;
+    dll_subroutine_update update_servo_with_limits ;
+    arraysizes_init  10 1 ;
+    arraysizes_update  5 9 ;
+	begin init ;
+      constant  1   3      ; Number of blades [-]
+      constant  2   1.0    ; Frequency of 2nd order servo model of pitch system [Hz]
+      constant  3   0.7    ; Damping ratio 2nd order servo model of pitch system [-]
+      constant  4  10.0    ; Max. pitch speed [deg/s]
+      constant  5  15.0    ; Max. pitch acceleration [deg/s^2]
+      constant  6  -5.0    ; Min. pitch angle [deg] 
+      constant  7  90.0    ; Max. pitch angle [deg] 	  
+	  constant  8 [Time pitch runaway]   ; Time for pitch runaway [s]
+	  constant  9 [Time stuck DLC22b]   ; Time for stuck blade 1 [s]
+	  constant 10 [Pitch 1 DLC22b]   ; Angle of stuck blade 1 [deg]
+	end init ;
+    begin output;
+      general time        ;  Time                         [s]     
+       dll inpvec 1 2     ;  Pitch1 demand angle          [rad]
+       dll inpvec 1 3     ;  Pitch2 demand angle          [rad]
+       dll inpvec 1 4     ;  Pitch3 demand angle          [rad]
+       dll inpvec 1 26    ;  Flag for emergency pitch stop         [0=off/1=on]
+    end output;           
+;
+    begin actions;    
+      constraint bearing2 angle pitch1 ; Angle pitch1 bearing    [rad]
+      constraint bearing2 angle pitch2 ; Angle pitch2 bearing    [rad]
+      constraint bearing2 angle pitch3 ; Angle pitch3 bearing    [rad]
+    end actions;                      
+  end type2_dll;
+;
+;	--- DLL for tower-blade tip distance -- ;
+  begin type2_dll; 
+    name towerclearance_mblade ; tower clearance DLL
+    filename  .\control\towerclearance_mblade ;
+    dll_subroutine_init initialize ;
+    dll_subroutine_update update ;
+    arraysizes_init  1 3 ;
+    arraysizes_update  12 6 ;
+    begin init ;    Variables passed into initialization function
+      constant  1 4.15  ; Tower radius at tower bottom [m] 
+      constant  2 2.75  ; Tower radius at tower top [m]
+      constant  3   3  ; Number of blade points to check [-]
+    end init ;
+    begin output;   Variables passed into update function
+      mbdy state pos  tower   1  0.00 global ; [1,2,3] global coordinates of tower base
+      mbdy state pos  tower  10  1.00 global ; [4,5,6] global coordinates of tower base
+      mbdy state pos blade1  26  1.00 global ; [7,8,9] global coordinates of point 1 (blade 1 100.0% R)
+      mbdy state pos blade2  26  1.00 global ; [10,11,12] global coordinates of point 2 (blade 2 100.0% R)
+      mbdy state pos blade3  26  1.00 global ; [12,13,14] global coordinates of point 3 (blade 3 100.0% R)
+    end output;    
+  end type2_dll;
+end dll;
+;---------------------------------------------------------------------------------------------------------------------------------------------------------------- 
+;
+begin output;
+  filename ./res/[Case folder]/[Case id.] ;
+  time [t0] [time stop] ;
+  data_format  [out_format];
+  buffer 1 ;
+;
+  general time;  
+  constraint bearing1 shaft_rot 2; angle and angle velocity 
+  constraint bearing2 pitch1 5;    angle and angle velocity 
+  constraint bearing2 pitch2 5;    angle and angle velocity 
+  constraint bearing2 pitch3 5;    angle and angle velocity 
+  aero omega ;
+  aero torque;
+  aero power;
+  aero thrust;
+  wind free_wind 1 0.0 0.0 -208.15; local wind at fixed position: coo (1=global,2=non-rotation rotor coo.), pos x, pos y, pos z
+  wind free_wind 1 0.0 0.0 -169; local wind at fixed position: coo (1=global,2=non-rotation rotor coo.), pos x, pos y, pos z
+  wind free_wind 1 0.0 0.0 -119; local wind at fixed position: coo (1=global,2=non-rotation rotor coo.), pos x, pos y, pos z
+  wind free_wind 1 0.0 0.0 -69; local wind at fixed position: coo (1=global,2=non-rotation rotor coo.), pos x, pos y, pos z  
+  wind free_wind 1 0.0 0.0 -29.85; local wind at fixed position: coo (1=global,2=non-rotation rotor coo.), pos x, pos y, pos z
+  ; Moments:
+  mbdy momentvec tower  1 1  tower # tower base ;
+  mbdy momentvec tower  10 2 tower # tower yaw bearing ;
+  mbdy momentvec shaft  4 1  shaft # main bearing ;
+  mbdy momentvec blade1 2  2 blade1 # blade 1 root ;
+  mbdy momentvec blade2 2  2 blade2 # blade 2 root ;
+  mbdy momentvec blade3 2  2 blade3 # blade 3 root ;
+  mbdy momentvec blade1 13 1 local # blade 1 50% local e coo ;
+  mbdy momentvec blade2 13 1 local # blade 2 50% local e coo ;
+  mbdy momentvec blade3 13 1 local # blade 3 50% local e coo ;
+  ; Displacements and accellerations
+  mbdy state pos tower 10 1.0 global only 1 # Tower top FA displ; 
+  mbdy state pos tower 10 1.0 global only 2 # Tower top SS displ; 
+  mbdy state acc tower 10 1.0 global only 1 # Tower top FA acc; 
+  mbdy state acc tower 10 1.0 global only 2 # Tower top SS acc; 
+;  
+  mbdy state pos blade1  26 1.0 blade1 # blade 1 tip pos ;
+  mbdy state pos blade2  26 1.0 blade2 # blade 2 tip pos ;
+  mbdy state pos blade3  26 1.0 blade3 # blade 3 tip pos ;
+  mbdy state pos blade1  26 1.0 global # gl blade 1 tip pos ;  
+; - Monitor Aerodynamics - ;
+  aero windspeed 3 1 1 72.5;
+  aero alfa 1 72.5; 
+  aero alfa 2 72.5; 
+  aero alfa 3 72.5; 
+  aero cl 1 72.5;
+  aero cl 2 72.5;
+  aero cl 3 72.5;
+  aero cd 1 72.5;
+  aero cd 2 72.5;
+  aero cd 3 72.5;    
+; - Main Controller - 
+; Output to controller
+  ; dll outvec 1 1 # time;
+  ; dll outvec 1 2 # slow speed shaft rad/s;
+  ; dll outvec 1 3 # pitch angle 1;
+  ; dll outvec 1 4 # pitch angle 2;
+  ; dll outvec 1 5 # pitch angle 3;
+  ; dll outvec 1 6 # WSP_x_global;
+  ; dll outvec 1 7 # WSP_y_global;
+  ; dll outvec 1 8 # WSP_z_global;
+  ; dll outvec 1 9 # Elec. pwr ;
+  ; dll outvec 1 10 # Grid flag ;
+; Input from controller
+  dll inpvec 1  1 # Generator torque reference            [Nm]   ;
+  dll inpvec 1  2 # Pitch angle reference of blade 1      [rad]  ;
+  dll inpvec 1  3 # Pitch angle reference of blade 2      [rad]  ;
+  dll inpvec 1  4 # Pitch angle reference of blade 3      [rad]  ;
+  ; dll inpvec 1  5 # Power reference                       [W]    ;
+  ; dll inpvec 1  6 # Filtered wind speed                   [m/s]  ;
+  ; dll inpvec 1  7 # Filtered rotor speed                  [rad/s];
+  ; dll inpvec 1  8 # Filtered rotor speed error for torque [rad/s];
+  ; dll inpvec 1  9 # Bandpass filtered rotor speed         [rad/s];
+  ; dll inpvec 1 10 # Proportional term of torque contr.    [Nm]   ;
+  ; dll inpvec 1 11 # Integral term of torque controller    [Nm]   ;
+  ; dll inpvec 1 12 # Minimum limit of torque               [Nm]   ;
+  ; dll inpvec 1 13 # Maximum limit of torque               [Nm]   ;
+  dll inpvec 1 14 # Torque limit switch based on pitch    [-]    ;
+  ; dll inpvec 1 15 # Filtered rotor speed error for pitch  [rad/s];
+  ; dll inpvec 1 16 # Power error for pitch                 [W]    ;
+  ; dll inpvec 1 17 # Proportional term of pitch controller [rad]  ;
+  ; dll inpvec 1 18 # Integral term of pitch controller     [rad]  ;
+  ; dll inpvec 1 19 # Minimum limit of pitch                [rad]  ;
+  ; dll inpvec 1 20 # Maximum limit of pitch                [rad]  ;
+  dll inpvec 1 21 # Torque reference from DT dammper      [Nm]  ;
+  dll inpvec 1 22 # Status signal                         [-]  ;
+  ; dll inpvec 1 23 # Total added pitch rate                [rad/s]  ;  
+  dll inpvec 1 24 # Filtered Mean pitch for gain sch      [rad]  ;  
+  dll inpvec 1 25 # Flag for mechnical brake              [0=off/1=on] ;
+  dll inpvec 1 26 # Flag for emergency pitch stop         [0=off/1=on] ;
+  dll inpvec 1 27 # LP filtered acceleration level        [m/s^2] ;
+; ; Output to generator model
+   ; dll outvec 2 1  # time ;
+   ; dll outvec 2 2  # Electrical torque reference [Nm] ;
+   ; dll outvec 2 3  # omega LSS ;
+; Input from generator model
+   dll inpvec 2 1  # Mgen LSS [Nm];
+   dll inpvec 2 2  # Pelec W ;
+   dll inpvec 2 3  # Mframe ;
+   dll inpvec 2 4  # Mgen HSS ;
+   dll inpvec 2 5  # Generator Pmech kW ;
+   dll inpvec 2 6  # Filtered Gen speed ;
+   dll inpvec 2 7  # Elec. pwr ;
+   dll inpvec 2 8  # Grid flag ;
+; Output to mechanical brake
+   dll inpvec 3 1 # Brake torque [Nm] ;
+; ; Input from mechanical brake
+   ; dll outvec 3 1 # Time [s] ;
+   ; dll outvec 3 2 # Generator LSS speed [rad/s] ;
+   ; dll outvec 3 3 # Deploy brake ;
+; ; Output to pitch servo
+   ; dll outvec 4 1 # time;
+   ; dll outvec 4 2 # pitchref 1;
+   ; dll outvec 4 3 # pitchref 2;
+   ; dll outvec 4 4 # pitchref 3;
+   ; dll outvec 4 5 # Emerg. stop;
+; Input from pitch servo
+   dll inpvec 4 1 # pitch 1;
+   dll inpvec 4 2 # pitch 2;
+   dll inpvec 4 3 # pitch 3;
+; Outputs from tower clearence DLL
+   dll inpvec 5 1 # Min clearance [m];
+   dll inpvec 5 2 # Idx closest pt to towr [-];
+   dll inpvec 5 3 # Rel dist tower bottom to top [-];
+   dll inpvec 5 4 # Tower rad at pt height [m];
+   dll inpvec 5 5 # Tower center x at pt height [m];
+   dll inpvec 5 6 # Tower center y at pt height [m];
+end output;
+;
+exit;
\ No newline at end of file
--- a/docs/update-conda-cluster.md
+++ b/docs/update-conda-cluster.md
-# Update conda ```wetb_py3``` environment and ```wetb```
+# Update conda ```py36-wetb``` environment and ```wetb```

 There are pre-configured miniconda/anaconda python environments installed on
 Gorm and Jess at:

 ```
-/home/python/miniconda3/envs/wetb_py3
+/home/python/miniconda3/envs/py36-wetb
 ```
 Note that these refer to the home drives of Gorm and Jess respectively and thus
 refer to two different directories (but are named the same).
@@ -15,20 +15,20 @@ Update the root Anaconda environment:
 conda update --all
 ```

-Activate the ```wetb_py3``` environment:
+Activate the ```py36-wetb``` environment:

 ```
-source activate wetb_py3
+source activate py36-wetb
 ```


-Update the ```wetb_py3``` environment:
+Update the ```py36-wetb``` environment:

 ```
 conda update --all
 ```

-Pull latest wetb changes and create re-distributable binary wheel package for ```wetb_py3```:
+Pull latest wetb changes and create re-distributable binary wheel package for ```py36-wetb```:

 ```
 cd /home/MET/repositories/tooblox/WindEnergyToolbox
@@ -45,5 +45,5 @@ pip install --no-deps -U dist/wetb-X.Y.Z.post0.devXXXXXXXX-cp35m-linux_x86_64.wh
 The option ```--no-deps``` is used here to avoid pip installing possible newer
 versions of packages that should be managed by conda. This only works when all
 dependencies of ```wetb``` are met (which is assumed by default for the
-```wetb_py3``` environment).
+```py36-wetb``` environment).

--- a/git_utils.py
+++ b/git_utils.py
+'''
+Created on 28. jul. 2017
+
+@author: mmpe
+'''
+import os
+import subprocess
+
+
+def _run_git_cmd(cmd, git_repo_path=None):
+    git_repo_path = git_repo_path or os.getcwd()
+    if not os.path.isdir(os.path.join(git_repo_path, ".git")):
+        raise Warning("'%s' does not appear to be a Git repository." % git_repo_path)
+    try:
+        process = subprocess.Popen(cmd,
+                                   stdout=subprocess.PIPE,
+                                   stderr=subprocess.PIPE,
+                                   universal_newlines=True,
+                                   cwd=os.path.abspath(git_repo_path))
+        stdout,stderr = process.communicate()
+        if process.returncode != 0:
+            raise EnvironmentError("%s\n%s"%(stdout, stderr))
+        return stdout.strip()
+
+    except EnvironmentError as e:
+        raise e
+        raise Warning("unable to run git")
+
+
+def get_git_branch(git_repo_path=None):
+    cmd = ["git", "rev-parse", "--abbrev-ref", "HEAD"]
+    return _run_git_cmd(cmd, git_repo_path)
+
+
+def get_git_version(git_repo_path=None):
+    cmd = ["git", "describe", "--tags", "--dirty", "--always"]
+    # format it will return: 'v0.1.0-12-g22668f0'
+    v = _run_git_cmd(cmd, git_repo_path)
+    # convert to something Pypi likes: 0.1.2.dev3.123456
+    # see also https://setuptools.pypa.io/en/latest/userguide/distribution.html
+    # and/or PEP440 https://peps.python.org/pep-0440/
+    v = v.replace('-', '.post', 1)
+    return v
+
+
+def get_tag(git_repo_path=None, verbose=False):
+    tag = _run_git_cmd(['git', 'describe', '--tags', '--always', '--abbrev=0'],
+                       git_repo_path)
+    if verbose:
+        print(tag)
+    return tag
+
+
+def set_tag(tag, push, git_repo_path=None):
+    _run_git_cmd(["git", "tag", tag], git_repo_path)
+    if push:
+        _run_git_cmd(["git", "push"], git_repo_path)
+        _run_git_cmd(["git", "push", "--tags"], git_repo_path)
+
+
+def update_git_version(version_module, git_repo_path=None):
+    """Update <version_module>.__version__ to git version"""
+
+    version_str = get_git_version(git_repo_path)
+    assert os.path.isfile(version_module.__file__)
+    with open(version_module.__file__, "w") as fid:
+        fid.write("__version__ = '%s'" % version_str)
+
+    # ensure file is written, closed and ready
+    with open(version_module.__file__) as fid:
+        fid.read()
+    return version_str
+
+
+def write_vers(vers_file='wetb/__init__.py', repo=None, skip_chars=1):
+    """Writes out version string as follows:
+        "last tag"-("nr commits since tag")-("branch name")-("hash commit")
+    and where nr of commits since last tag is only included if >0,
+    branch name is only inlcuded when not on master,
+    and hash commit is only included when not at a tag (when nr of commits > 0)
+    """
+    if not repo:
+        repo = os.getcwd()
+    version_long = get_git_version(repo)
+    branch = get_git_branch(repo)
+
+    verel = version_long.split('-')
+    # tag name
+    version = verel[0][skip_chars:]
+    # number of commits since last tag, only if >0
+    nr_commits = 0
+    if len(verel) > 1:
+        try:
+            nr_commits = int(verel[1])
+        except ValueError:
+            nr_commits = -1
+        if nr_commits > 0:
+            version += '-' + verel[1]
+    # branch name, only when NOT on master
+    if branch != 'master':
+        version += '-' + branch
+    # hash commit, only if not at tag
+    if len(verel) > 2 and nr_commits > 0:
+        # first character on the hash is always a g (not part of the hash)
+        version += '-' + verel[2][1:]
+    # if "-HEAD" is added to the version, which pypi does not like:
+    if version.endswith('-HEAD'):
+        version = version[:-5]
+    print(version_long)
+    print('Writing version: {} in {}'.format(version, vers_file))
+
+    with open(vers_file, 'r') as f:
+        lines = f.readlines()
+    for n, l in enumerate(lines):
+        if l.startswith('__version__'):
+            lines[n] = "__version__ = '{}'\n".format(version)
+    for n, l in enumerate(lines):
+        if l.startswith('__release__'):
+            lines[n] = "__release__ = '{}'\n".format(version)
+    with open(vers_file, 'w') as f:
+        f.write(''.join(lines))
+    return version
+
+
+def rename_dist_file():
+    for f in os.listdir('dist'):
+        if f.endswith('whl'):
+            split = f.split('linux')
+            new_name = 'manylinux1'.join(split)
+            old_path = os.path.join('dist', f)
+            new_path = os.path.join('dist', new_name)
+            os.rename(old_path, new_path)
+
+
+def main():
+    """Example of how to run (pytest-friendly)"""
+    if __name__ == '__main__':
+        pass
+
+
+main()
--- a/overview.xlsx
+++ b/overview.xlsx
--- a/pyproject.toml
+++ b/pyproject.toml
+[build-system]
+requires = [
+    "setuptools>=60",
+    "setuptools-scm>=8.0"]
+build-backend = "setuptools.build_meta"
+
+[project]
+name = "wetb"
+authors = [{name="DTU Wind and Energy Systems"}]
+description = "The Wind Energy Toolbox (or ```wetb```, pronounce as wee-tee-bee) is a collection of Python scripts that facilitate working with (potentially a lot) of HAWC2, HAWCStab2, FAST or other text input based simulation tools."
+dependencies = [
+    'certifi',
+    'click',
+    'Cython',
+    'h5py',
+    'Jinja2',
+    'lxml',
+    'matplotlib',
+    'pillow',
+    'mock',
+    'numpy',
+    'numba',
+    'openpyxl',
+    'pandas',
+    'paramiko',
+    'psutil',
+    'pytest',
+    'pytest-cov',
+    'scipy',
+    'sshtunnel',
+    'tables',
+    'tqdm',
+    'xarray',
+    'xlwt',
+    'XlsxWriter',
+]
+license = {text = "wetb is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License (GPL, http://www.gnu.org/copyleft/gpl.html) as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.  wetb is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more details http://www.gnu.org/licenses/  We encourage you to submit new code for possible inclusion in future versions of wetb."}
+dynamic = ["version"]
+
+[project.urls]
+repository = "https://gitlab.windenergy.dtu.dk/toolbox/WindEnergyToolbox"
+documentation = "https://toolbox.pages.windenergy.dtu.dk/WindEnergyToolbox/"
+
+
+[project.optional-dependencies]
+prepost = ["openpyxl", "tables", "xlwt", "Cython"]
+all = ["openpyxl", "tables", "xlwt", "Cython", "paramiko", "sshtunnel", 'pytest', 'mock', 'click']
+
+[tool.setuptools_scm]
+version_scheme = "no-guess-dev"
+
+[tool.setuptools]
+packages = ["wetb"]
--- a/remove-py2-support.sh
+++ b/remove-py2-support.sh
+#!/bin/bash
+
+find ./ -type f -iname "*.py" -exec sed -i '/from __future__ import print_function/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from __future__ import division/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from __future__ import unicode_literals/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from __future__ import absolute_import/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from io import open/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import object/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import map/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import chr/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import dict/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import super/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import zip/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import range/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import str/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from builtins import int/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/from future import standard_library/d' {} \;
+find ./ -type f -iname "*.py" -exec sed -i '/standard_library.install_aliases()/d' {} \;
+
--- a/requirements.txt
+++ b/requirements.txt
-# Add your requirements here like:
-six
-numpy>=1.4
-scipy>=0.9
-matplotlib
-pytest
-xlrd
-xlwt
-openpyxl
-h5py
-pandas
-tables
-future
-paramiko
-psutil
-pbr
-PyScaffold
-pytest-cov
-
--- a/setup.cfg
+++ b/setup.cfg
-[metadata]
-name = wetb
-summary = Wind Energy Toolbox
-author = DTU Wind Energy
-author-email = mmpe@dtu.dk
-license = GPLv3
-home-page = https://gitlab.windenergy.dtu.dk/toolbox/WindEnergyToolbox
-description-file = README
-# Add here all kinds of additional classifiers as defined under
-# https://pypi.python.org/pypi?%3Aaction=list_classifiers
-classifiers = Development Status :: 4 - Beta,
-              Programming Language :: Python,
-              Programming Language :: Python :: 2.7,
-              Programming Language :: Python :: 3,
-              Programming Language :: Python :: 3.3,
-              Programming Language :: Python :: 3.4,
-			  Programming Language :: Python :: 3.5,
-              Environment :: Console,
-              Intended Audience :: Education,
-              Intended Audience :: Science/Research,
-              License :: OSI Approved :: GPL License,
-              Operating System :: OS Independent,
-              Operating System :: POSIX :: Linux,
-              Operating System :: Unix,
-              Operating System :: MacOS,
-              Operating System :: Microsoft :: Windows
-              Topic :: Scientific/Engineering :: Mathematics
-
-[entry_points]
-# Add here console scripts like:
-# console_scripts =
-#     hello_world = wetb.module:function
-# as well as other entry_points.
-
-
-[files]
-# Add here 'data_files', 'packages' or 'namespace_packages'.
-# Additional data files are defined as key value pairs of source and target:
-packages =
-    wetb
-# data_files =
-#    share/wetb_docs = docs/*
-
-[extras]
-# Add here additional requirements for extra features, like:
-# PDF =
-#    ReportLab>=1.2
-#    RXP
-#ALL =
-#    django
-#    cookiecutter
-
-[test]
-# py.test options when running `python setup.py test`
-#addopts = tests
-
-[pytest]
-# Options for py.test:
-# Specify command line options as you would do when invoking py.test directly.
-# e.g. --cov-report html (or xml) for html/xml output or --junitxml junit.xml
-# in order to write a coverage file that can be read by Jenkins.
-#addopts =
-#    --cov wetb --cov-report term-missing
-#    --verbose
-python_files = WindEnergyToolbox/wetb/*
-
-[aliases]
-docs = build_sphinx
-
-[bdist_wheel]
-# Use this option if your package is pure-python
-universal = 0
-
-[build_sphinx]
-# Options for Sphinx build
-source_dir = docs
-build_dir = docs/_build
-
-[pbr]
-# Let pbr run sphinx-apidoc
-autodoc_tree_index_modules = True
-# autodoc_tree_excludes = ...
-# Let pbr itself generate the apidoc
-# autodoc_index_modules = True
-# autodoc_exclude_modules = ...
-# Convert warnings to errors
-# warnerrors = True
-
-[devpi:upload]
-# Options for the devpi: PyPI serer and packaging tool
-# VCS export must be deactivated since we are using setuptools-scm
-no-vcs = 1
-formats = bdist_wheel
--- a/setup.py
+++ b/setup.py
-#!/usr/bin/env python
-# -*- coding: utf-8 -*-
-"""
-    Setup file for wetb.
-
-    This file was generated with PyScaffold 2.5, a tool that easily
-    puts up a scaffold for your new Python project. Learn more under:
-    http://pyscaffold.readthedocs.org/
-"""
-
-import os
-import sys
-from setuptools import setup
-
-try:
-    from pypandoc import convert_file
-    read_md = lambda f: convert_file(f, 'rst')
-except ImportError:
-    print("warning: pypandoc module not found, could not convert Markdown to RST")
-    read_md = lambda f: open(f, 'r').read()
-
-import numpy as np
-from distutils.extension import Extension
-from Cython.Distutils import build_ext
-
-
-def setup_package():
-
-    path = 'wetb/fatigue_tools/rainflowcounting/'
-    module = 'wetb.fatigue_tools.rainflowcounting'
-    names = ['pair_range', 'peak_trough', 'rainflowcount_astm']
-    extlist = [Extension('%s.%s' % (module, n),
-                         [os.path.join(path, n)+'.pyx'],
-                         include_dirs=[np.get_include()]) for n in names]
-
-    needs_sphinx = {'build_sphinx', 'upload_docs'}.intersection(sys.argv)
-    sphinx = ['sphinx'] if needs_sphinx else []
-    setup(setup_requires=['six', 'pyscaffold>=2.5a0,<2.6a0'] + sphinx,
-          cmdclass = {'build_ext': build_ext},
-          ext_modules = extlist,
-          use_pyscaffold=True,
-          long_description=read_md('README.md'))
-
-
-if __name__ == "__main__":
-    setup_package()
--- a/tests/__init__.py
+++ b/tests/__init__.py
+import numpy as np
+npt = np.testing
\ No newline at end of file
--- a/tests/run_main.py
+++ b/tests/run_main.py
+import os
+import sys
+
+from unittest import mock
+import pytest
+import matplotlib.pyplot as plt
+
+
+def run_module_main(module):
+    # check that all main module examples run without errors
+    if os.name == 'posix' and "DISPLAY" not in os.environ:
+        pytest.xfail("No display")
+
+    def no_show(*args, **kwargs):
+        pass
+    plt.show = no_show  # disable plt show that requires the user to close the plot
+
+    def no_print(s):
+        pass
+
+    try:
+        with mock.patch.object(module, "__name__", "__main__"):
+            with mock.patch.object(module, "print", no_print):
+                getattr(module, 'main')()
+    except Exception as e:
+        raise type(e)(str(e) +
+                      ' in %s.main' % module.__name__).with_traceback(sys.exc_info()[2])
--- a/tests/test_examples.py
+++ b/tests/test_examples.py
+import importlib
+import os
+import pkgutil
+import warnings
+import mock
+import pytest
+import matplotlib.pyplot as plt
+import sys
+from wetb import examples
+from tests.run_main import run_module_main
+
+
+def get_main_modules():
+    package = examples
+    modules = []
+    for _, modname, _ in pkgutil.walk_packages(package.__path__, package.__name__ + '.'):
+        with warnings.catch_warnings():
+            warnings.simplefilter("ignore")
+            m = importlib.import_module(modname)
+
+        if 'main' in dir(m):
+            modules.append(m)
+    return modules
+
+
+def print_main_modules():
+    print("\n".join([m.__name__ for m in get_main_modules()]))
+
+
+@pytest.mark.parametrize("module", get_main_modules())
+def test_main(module):
+    run_module_main(module)
+
+
+if __name__ == '__main__':
+    print_main_modules()
--- a/wetb/__init__.py
+++ b/wetb/__init__.py
-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import division
-from __future__ import absolute_import
-from future import standard_library
-standard_library.install_aliases()
 test = "TEST"
-try:
-    import pkg_resources
-    __version__ = pkg_resources.safe_version(pkg_resources.get_distribution(__name__).version)
-except:
-    __version__ = 'unknown'
+__version__ = '0.0.10'
+# try:
+#     import pkg_resources
+#     __version__ = pkg_resources.safe_version(pkg_resources.get_distribution(__name__).version)
+# except:
+#     __version__ = 'unknown'
--- a/wetb/bladed/prj2hawc.py
+++ b/wetb/bladed/prj2hawc.py
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Wed Jun 23 11:20:28 2021
+
+@author: dave
+"""
+
+import os
+from os.path import join as pjoin
+
+import numpy as np
+import pandas as pd
+
+# from wetb.prepost import misc
+from wetb.hawc2 import (HTCFile, AEFile, PCFile, StFile)
+
+from wetb.bladed.readprj import ReadBladedProject
+
+
+# TODO: move to misc?
+def get_circular_inertia(d1, t, m):
+    d2 = d1 - 2*t
+    A = np.pi/4*(d1**2-d2**2)
+    # area moment of inertia x=y
+    I_g = np.pi/64*(d1**4-d2**4)
+    # mass moment of inertia
+    I_m = I_g/A*m
+    return A, I_g, I_m
+
+
+# TODO: move to HTCFile?
+def add_c2_arr_htc(htc, c2_arr, body):
+    """add a c2_def in the form of an array to a htc c2_def section. Needs to
+    be an existing body, and the current c2_def will be overwritten.
+    """
+    nr_nodes = c2_arr.shape[0]
+    c2_def = htc['new_htc_structure'].get_subsection_by_name(body)['c2_def']
+
+    keys = ['sec'] + [f'sec__{i}' for i in range(2, nr_nodes+1)]
+
+    c2_def['nsec'].values[0] = nr_nodes
+    for i, key in enumerate(keys):
+        # control the formatting, otherwise it looks very messy in the htc file
+        # this we do by converting the values to strings
+        vals_fmt = [f'{k: 15.07f}' for k in c2_arr[i,:].tolist()]
+        # node nr is stripped from leading spaces, so it doesn't help adding them
+        values = [i+1] + vals_fmt
+        if key in c2_def:
+            c2_def[key].values[:] = values[:]
+        else:
+            c2_def.add_line(name=key, values=values)
+
+    return htc
+
+
+# TODO: move to HTCFile?
+def get_c2def(htc, body_name):
+    """
+    """
+
+    # select the c2def section from the shaft
+    c2_def = htc['new_htc_structure'].get_subsection_by_name(body_name)['c2_def']
+    # safe all c2_def positions in one array
+    c2_arr = np.zeros((len(c2_def.keys())-1,4))
+    # iterate over all nsec lines, but ignore the line that gives the number of points
+    for i, sec in enumerate(c2_def.keys()[1:]):
+        c2_arr[i,:] = c2_def[sec].values[1:]
+
+    return pd.DataFrame(c2_arr, columns=['x', 'y', 'z', 'twist'])
+
+
+class Convert2Hawc(ReadBladedProject):
+    """Based on the BLADED user manual v4.8
+    """
+
+    def __init__(self, fname):
+        super().__init__(fname)
+
+    def get_pc(self, setname_filter=None):
+        # ---------------------------------------------------------------------
+        # we ignore the airfoils blocks at the start of the ADAT sub-section
+        # ---------------------------------------------------------------------
+        # nr_airfoils = prj.get_key('ADAT', 'NFOILS')[0,0]
+        # airfoil_keys = ['NFOIL', 'FTYPE', 'NTHICK', 'NREYN']
+        # airfoils = {}
+        # for k in range(nr_airfoils):
+        #     key_app = ''
+        #     if k > 0:
+        #         key_app = f'__{k+1}'
+        #     airfoils[f'NFOIL_{k+1}'] = {}
+        #     for key in airfoil_keys:
+        # ---------------------------------------------------------------------
+
+        # total number of pc sets
+        nsets = self.get_key('ADAT', 'NSETS')[0,0]
+        key_app = [''] + [f'__{k}' for k in range(2, nsets+1)]
+
+        # FIXME: are all sets always named consistantly? Meaning only the
+        # thickness changes in the setname, and the rest remains the same?
+        # as a work-around: ignore sets that contains setname_filter in SETNAME
+        if setname_filter is not None:
+            drop = []
+            for i, ka in enumerate(key_app):
+                setname = self.get_key('ADAT', f'SETNAME{ka}')[0,0]
+                if setname.find(setname_filter) > -1:
+                    drop.append(ka)
+            key_app = set(key_app) - set(drop)
+            # update nsets
+            nsets = len(key_app)
+
+        # initiate PCFile objects
+        pc = PCFile()
+        tc_arr = np.ndarray((nsets), dtype=np.float32)
+        pc_lst = []
+
+        for i, ka in enumerate(key_app):
+
+            # and assure cm is around quarter chord point
+            assert self.get_key('ADAT', 'XA')[0,0]==25
+
+            # prj.get_key('ADAT', f'SETNAME{ka}')
+            # prj.get_key('ADAT', f'SETNB{ka}')
+            # prj.get_key('ADAT', f'COMMENTS{ka}')
+            # prj.get_key('ADAT', f'DATE{ka}')
+            # prj.get_key('ADAT', f'REYN{ka}')
+            # prj.get_key('ADAT', f'DEPL{ka}')
+            # prj.get_key('ADAT', f'NALPHA{ka}')
+            # prj.get_key('ADAT', f'NVALS{ka}')
+
+            tc_arr[i] = self.get_key('ADAT', f'THICK{ka}')
+            shape = (self.get_key('ADAT', f'NALPHA{ka}')[0,0], 4)
+            pcs = np.ndarray(shape, dtype=np.float32)
+            pcs[:,0] = self.get_key('ADAT', f'ALPHA{ka}')[0,:]
+            pcs[:,1] = self.get_key('ADAT', f'CL{ka}')[0,:]
+            pcs[:,2] = self.get_key('ADAT', f'CD{ka}')[0,:]
+            pcs[:,3] = self.get_key('ADAT', f'CM{ka}')[0,:]
+            pc_lst.append(pcs)
+
+        # in BGEOM the airfoil number is referred to and not the thickness
+        # sort properly on thickness
+        isort = tc_arr.argsort()
+        pc.pc_sets = {1:(tc_arr[isort], [pc_lst[k] for k in isort])}
+
+        return pc
+
+    def get_ae(self):
+        # BGEOMMB defines for each element the start and end point. We only
+        # need half that of course.
+        chord = self.get_key('BGEOMMB', 'CHORD')[0,0::2]
+        # FIXME: is it radius or curved lenght
+        radius = self.get_key('BGEOMMB', 'RJ')[0,0::2]
+        thickness = self.get_key('BGEOMMB', 'BTHICK')[0,0::2]
+        # FIXME: multiple pc sets can be defined but are ignored for now.
+        pc_set_id = np.ones(chord.shape)
+
+        ae = AEFile()
+        ae.add_set(radius, chord, thickness, pc_set_id, set_id=None)
+
+        return ae
+
+    def get_blade_c2_def(self):
+
+        # BLADENAME
+        # DISTMODE	RJ
+        # NBE	100
+        # FOIL
+        # MOVING
+
+        # CE_X
+        # CE_Y
+        diam = self.get_key('RCON', 'DIAM')[0,0]
+        cone = self.get_key('RCON', 'CONE')[0,0]
+        tilt = self.get_key('RCON', 'TILT')[0,0]
+
+        # in BLADED all is wrt LE, so we need chord
+        chord = self.get_key('BGEOMMB', 'CHORD')[0,0::2]
+
+        # c2_def = htc['new_htc_structure'].get_subsection_by_name('blade1')['c2_def']
+        twist = self.get_key('BGEOMMB', 'TWIST')[0,0::2]
+
+        # FIXME: what is the difference between RJ and DIST?
+        # radius can be given eiter as curved length or radius
+        # Distance: This can be entered as a distance along the blade or as a
+        # distance along the blade root Z-axis. Select the appropriate option
+        # at the base of the screen.
+        dist = self.get_key('BGEOMMB', 'DIST')[0,0::2]
+        rj = self.get_key('BGEOMMB', 'RJ')[0,0::2]
+        # reference frame is EC, that seems to be indicated on the plots of the
+        # blade rotor coordinates/pitch/axis tc (fig 4-1 and 4-2)
+        ref_x = self.get_key('BGEOMMB', 'REF_X')[0,0::2]
+        ref_y = self.get_key('BGEOMMB', 'REF_Y')[0,0::2]
+        # EC wrt to LE in percentage of chord
+        ce_x = self.get_key('BGEOMMB', 'CE_X')[0,0::2]
+        ce_y = self.get_key('BGEOMMB', 'CE_Y')[0,0::2]
+
+        # HAWC2 variables
+        theta = twist*-1
+        # theta_d = theta*180/np.pi
+        ea_x = (50-ce_y)*0.01*chord
+        ea_y = ce_x*0.01*chord
+
+        # safe all c2_def positions in one array
+        c2_arr = np.zeros((len(twist),4))
+        c2_arr[:,0] = -1*(ref_y + ((ea_x*np.cos(theta) - ea_y*np.sin(theta))))
+        c2_arr[:,1] =    (ref_x - ((ea_x*np.sin(theta) + ea_y*np.cos(theta))))
+        # assume DIST is curved length
+        c2_arr[:,2] = rj #dist
+        c2_arr[:,3] = theta*180/np.pi
+
+        return c2_arr
+
+    def get_blade_st(self):
+        # BMASSMB
+        # BSTIFFMB
+
+        # in HAWC2 structural pitch is relative to the twist
+        twist = self.get_key('BGEOMMB', 'TWIST')[0,0::2]
+        dist = self.get_key('BGEOMMB', 'DIST')[0,0::2]
+        chord = self.get_key('BGEOMMB', 'CHORD')[0,0::2]
+        thickness = self.get_key('BGEOMMB', 'BTHICK')[0,0::2]
+
+        # CG location, in percentage chord wrt LE
+        cm_x = self.get_key('BMASSMB', 'CM_X')[0,0::2]
+        cm_y = self.get_key('BMASSMB', 'CM_Y')[0,0::2]
+        mass = self.get_key('BMASSMB', 'MASS')[0,0::2]
+        siner = self.get_key('BMASSMB', 'SINER')[0,0::2]
+
+        # Radii of gyration ratio: the radius of gyration of mass about y_m
+        # divided by the radius of gyration of mass about x_m. This defaults to
+        # the relative profile thickness but can be over-written by un-checking
+        # the “Use default radii of gyration ratio” checkbox.
+        rgratio = self.get_key('BMASSMB', 'RGRATIO')[0,0::2]
+        # FIXME: DEF_RGRATIO but which value belongs to which option?
+        # assume 0 means "default is off"
+        def_rgratio = self.get_key('BMASSMB', 'DEF_RGRATIO')[0,0]
+
+        # FIXME: BETA_M: figure 3-2 from v4.8: mass axis orientation, radians?
+        beta_m = self.get_key('BMASSMB', 'BETA_M')[0,0::2]
+        # Mass axis orientation: the orientation of the principle axis of inertia.
+        # This defaults to the orientation of aerodynamic twist, but can be
+        # over-written by un-checking the “Use default mass axis orientation”
+        # checkbox. (See diagram below)
+        # FIXME: DEF_BETA_M but which value belongs to which option?
+        # assume 0 means "default is off"
+        def_beta_m = self.get_key('BMASSMB', 'DEF_BETA_M')[0,0]
+
+        # DEF_RGRATIO	0
+        # DEF_BETA_M	0
+
+        # RADPOS
+        # PITCHPOS
+
+        # DEF_RGRATIO	0
+        # DEF_BETA_M	0
+        # ICEDBLADES	 0, 0, 0
+        # ICEDENSITY	 700
+        # TIPCHORD	 0
+        # NPOINT	3
+        # RADPOS	 .2, 30.12345, 40.12345
+        # PITCHPOS	 .2, 30.1, 40.1
+        # PM_X	 0, 0, 0
+        # PM_Y	 20, 10, 20
+        # PMASS	 1000, 30, 1
+
+        # elastic axis, in percentage chord wrt LE
+        ce_x = self.get_key('BGEOMMB', 'CE_X')[0,0::2]
+        ce_y = self.get_key('BGEOMMB', 'CE_Y')[0,0::2]
+
+        # 3.5 Blade stiffness distribution (page 19)
+        # The stiffness must be defined about the principal axis of inertia at
+        # each blade station (see 3.1). The stiffness is the product of Young’s
+        # Modulus for the material and the second moment of area for the xp or
+        # yp directions as appropriate. The principal axis orientation is defined
+        # as an input, and defaults to the aerodynamic twist. In this case it is
+        # assumed to be parallel and perpendicular to the chord line. If the
+        # principal axis orientation is different from the aerodynamic twist, click
+        # the Use default principal axis orientation to off. (see diagram below)
+        eiflap = self.get_key('BSTIFFMB', 'EIFLAP')[0,0::2]
+        eiedge = self.get_key('BSTIFFMB', 'EIEDGE')[0,0::2]
+        beta_s = self.get_key('BSTIFFMB', 'BETA_S')[0,0::2]
+        gj = self.get_key('BSTIFFMB', 'GJ')[0,0::2]
+        cs_x = self.get_key('BSTIFFMB', 'CS_X')[0,0::2]
+        cs_y = self.get_key('BSTIFFMB', 'CS_Y')[0,0::2]
+        gaflap = self.get_key('BSTIFFMB', 'GAFLAP')[0,0::2]
+        gaedge = self.get_key('BSTIFFMB', 'GAEDGE')[0,0::2]
+        # FIXME: DEF_BETA_S probably says what the frame of reference is
+        # assume 0 means "default is off"
+        # self.get_key('BMASSMB', 'DEF_BETA_S')[0,0::2]
+
+        # for plotting purposes
+        sel = beta_m < -1
+        beta_m[sel] += np.pi
+        sel = beta_s < -1
+        beta_s[sel] += np.pi
+
+        # calculate the curved length
+        # c2_arr = self.get_blade_c2_def()
+        # curved_len_arr = curved_length(c2_arr)
+
+        ea_x = (50-ce_y)*0.01*chord
+        ea_y = ce_x*0.01*chord
+        cg_x = (50-cm_y)*0.01*chord
+        cg_y = cm_x*0.01*chord
+
+        # mass moments
+        Ipm_y = siner + (mass*(cg_x-ea_x)**2)
+        ri_y = np.sqrt( (Ipm_y / (1+rgratio**2)) / mass )
+        ri_x = rgratio * ri_y
+        # in BLADED rgratio = y/x, but in HAWC2 it x=y so we swap them
+
+        # just choose E, G so the rest follows
+        E, G = 1e10, 1e09
+        # however, we don't know EA (probably because stiff?), so assume A is
+        # a full box based on chord length and maximum airfoil thickness
+        # this will make the extensional stiffness EA quite large, but the user
+        # would still end up with a physical meaningful value
+        A = chord*chord*thickness/100
+
+        starr = np.ndarray((len(dist),19))
+        starr[:,0] = dist
+        starr[:,1] = mass
+        starr[:,2] = cg_x
+        starr[:,3] = cg_y
+        # radius of gyration
+        starr[:,4] = ri_x
+        starr[:,5] = ri_y
+        # shear center
+        starr[:,6] = (50-cs_y)*0.01*chord
+        starr[:,7] = cs_x*0.01*chord
+        # E, G: choose
+        starr[:,8] = E
+        starr[:,9] = G
+        starr[:,10] = eiflap/E
+        starr[:,11] = eiedge/E
+        starr[:,12] = gj/G
+        # shear factor
+        starr[:,13] = gaflap/(G*A)
+        starr[:,14] = gaedge/(G*A)
+        starr[:,15] = A
+        starr[:,16] = -1*(beta_s - twist) # structural pitch
+        starr[:,17] = ea_x
+        starr[:,18] = ea_y
+
+        st = StFile()
+        st.main_data_sets = {1:{1:starr}}
+        # st.cols = wetb.hawc2.st_file.stc.split()
+        return st
+
+    def get_tower_c2_def_st(self):
+        """Ignores nodes that are not at the centerline
+
+        Returns
+        -------
+        c2_arr : TYPE
+            DESCRIPTION.
+        st : TYPE
+            DESCRIPTION.
+
+        """
+
+        t_n_e = int(self.get_key('TGEOM', 'NTE'))
+        t_n_n = int(self.get_key('TGEOM', 'NTS'))
+        t_el_mat = self.get_key('TGEOM', 'ELSTNS')
+        t_n_coord = self.get_key('TGEOM', 'TCLOCAL')
+        # tdiam = self.get_key('TGEOM', 'SEAL')
+        # tdiam = self.get_key('TGEOM', 'HYDRO')
+        towmgt = self.get_key('TGEOM', 'TOWMGT')
+        # MATERIAL, NOMATS
+        material = self.get_key('TMASS', 'MATERIAL')[0,0::2]
+        nomats = self.get_key('TMASS', 'NOMATS')
+
+        # Index of tower nodes (nodes at zero x-y positions)
+        t_a = np.where((abs(t_n_coord[:,0])+abs(t_n_coord[:,1])) == 0)[0]
+        t_nodes = t_n_coord[t_a]
+        # sort index according to z position of the nodes
+        n_ind = t_a[np.argsort(t_nodes[:,2])]
+        t_nodes = t_n_coord[n_ind]
+        # Index for tower elemetns
+        # Elemet indexes which does not contain both tower nodes
+        ind_not = [i for i,j in enumerate(t_el_mat) if not (j[0]-1) in n_ind
+                   or not (j[1]-1) in n_ind]
+        # The elemetn indexes excluding the ind_not
+        ind_init = np.array([i for i in range(t_el_mat.shape[0]) if not i in ind_not])
+        sort_el = []
+        for i in n_ind[:-1]:
+            sort_el.append(np.where(t_el_mat[ind_init][:,0] == i+1)[0][0])
+        ind_el = ind_init[sort_el]
+
+        # safe all c2_def positions in one array
+        c2_arr = np.zeros((t_nodes.shape[0],4))
+        c2_arr[:,0:2] = t_nodes[:,:-1].copy()
+        c2_arr[:,2] = -t_nodes[:,-1] + t_nodes[0,-1] # It is minus for HAWC2 model
+        c2_arr[:,3] = 0.0
+        # element properties
+        t_d =  self.get_key('TGEOM', 'TDIAM')[ind_el]
+        material = self.get_key('TMASS', 'MATERIAL')[0,0::2][ind_el]
+        t_thick = self.get_key('TMASS', 'WALLTHICK')[ind_el]
+        towmgt = self.get_key('TGEOM', 'TOWMGT')[ind_el]
+        # Mass
+        t_m = self.get_key('TMASS', 'TOWM')[ind_el]
+        t_p_m = self.get_key('TMASS', 'TOWPMI')[ind_el] # point mass for elements and nodes
+        t_flood = self.get_key('TMASS', 'FLOODED')[ind_el]
+        # material props from prj file
+        name_mat = np.unique(material)
+        mat_prop = {}
+        for i in name_mat:
+            mat_prop[i] = self.get_key('TMASS', "%s"%i) # rho, E, G
+
+        t_mat = np.array([mat_prop[i][0] for i in material])
+
+        if np.sum(t_thick[1:,0] - t_thick[:-1,1]) > 1e-9:
+            print('Diameters are not continous along the tower')
+        if (name_mat.shape[0] > 1):
+            print('There are %i materials in tower!' %name_mat.shape[0])
+        d = np.zeros(t_d.shape[0]+1)
+        t = np.zeros(t_thick.shape[0]+1)
+        mat = np.zeros((t_mat.shape[0]+1,t_mat.shape[1]))
+        d[1:] = t_d[:,1]
+        t[1:] = t_thick[:,1]
+        mat[1:,:] = t_mat.copy()
+        d[0] = t_d[0,1]
+        t[0] = t_thick[0,1]
+        mat[0,:] = t_mat[0,:].copy()
+
+        A, I_g, I_m = get_circular_inertia(d, t, mat[:,0])
+        # FIXME: MAYBE ANOTHER FIX FOR DIFFERENT MATERIALS
+        # return c2_arr
+        starr = np.zeros((c2_arr.shape[0],19))
+        starr[:,0] = -c2_arr[:,2] #
+        starr[:,1] = A*mat[:,0]
+        starr[:,2] = 0.0 # no cg offset
+        starr[:,3] = 0.0 # no cg offset
+        # radius of gyration
+        starr[:,4] = np.sqrt(I_m/mat[:,0]) # ri_x = (I_m/m)**0.5 = (I_g/A)**0.5
+        starr[:,5] = np.sqrt(I_m/mat[:,0]) # ri_y = (I_m/m)**0.5 = (I_g/A)**0.5
+        # shear center
+        starr[:,6] = 0.0 # no shear center offset
+        starr[:,7] = 0.0 # no shear center offset
+        # E, G: choose
+        starr[:,8] = mat[:,1]
+        starr[:,9] = mat[:,2]
+        starr[:,10] = I_g
+        starr[:,11] = I_g
+        starr[:,12] = I_g*2
+        # shear factor
+        starr[:,13] = 3/4 # for circular section check it again
+        starr[:,14] = 3/4 # for circular section check it again
+        starr[:,15] = A
+        starr[:,16] = 0.0 # structural pitch
+        starr[:,17] = 0.0
+        starr[:,18] = 0.0
+
+        st = StFile()
+        st.main_data_sets = {1:{1:starr}}
+
+        return c2_arr, st
+
+    def get_hub(self):
+        # Blade root length, in blade root coordinates (with tilt and cone)
+        # see figure 4-1 in Blade User Manual 4.8
+        root = self.get_key('RCON', 'ROOT')[0,0]
+
+        # Hub mass: the mass of the hub, including the spinner and any blade
+        #   root section.
+        # Hub mass centre: the distance from the intersection of the shaft and
+        #   blade axes to the centre of mass of the hub, in a direction measured
+        #   away from the tower.
+        # Moments of inertia: the moment of inertia of the hub mass about the
+        #   shaft axis must be defined.
+        # The inertia about an axis perpendicular to the shaft may also be
+        #   entered with its origin about the hub centre of mass.
+
+        hubmas = self.get_key('RMASS', 'HUBMAS')[0,0]
+        hubine = self.get_key('RMASS', 'HUBINE')[0,0]
+        # FIXME: which axis is this? Would that translate to both xx and yy in H2?
+        hubine2 = self.get_key('RMASS', 'HUBINE2')[0,0]
+
+        # Enter the Spinner diameter. This is the diameter of any spinner or
+        # nose-cone, within which the blades themselves experience no
+        # aerodynamic forces.
+        spind = self.get_key('RCON', 'SPIND')[0,0]
+
+        cmass_hub = {'x':0, 'y':0, 'z':0, 'm':hubmas,
+                     'Ixx':0, 'Iyy':0, 'Izz':hubine}
+
+        return cmass_hub, root
+
+    def get_drivertain(self, len_shaft):
+
+        # TODO: brake
+        # MSTART BRAKE
+
+        # DTRAIN
+        # dtrain = self.get_key('DTRAIN')
+        ginert = self.get_key('DTRAIN', 'GINERT')[0,0] # generator inertia
+        gratio = self.get_key('DTRAIN', 'GRATIO')[0,0]
+        # The additional inertia of the high speed shaft may also be specified
+        # along with the inertia of the gearbox which is referred to the HSS.
+        gbxinert = self.get_key('DTRAIN', 'GBXINERT')[0,0] # gearbox inertia
+        # FIXME: brake position?
+        bpos = self.get_key('DTRAIN', 'BPOS')[0,0]
+        # LSS seems to have a DOF indicator in LSSDOF, and has 6 elements
+        # 4th element is torsion it seems
+        klss = self.get_key('DTRAIN', 'KLSS')[0,3]
+        # FIXME: what does this damping value mean? Simply C in the dyn system
+        # theta_dt_dt*I + C*theta_dt + K*theta = 0 ?
+        dlss = self.get_key('DTRAIN', 'DLSS')[0,3]
+        # HSS only has one DOF
+        khss = self.get_key('DTRAIN', 'KHSS')[0,0]
+        dhss = self.get_key('DTRAIN', 'DHSS')[0,0]
+
+        # Convert torsional stiffness definition into a beam stiffness
+        # that applies for the considered shaft length (c2_def)
+        G = 80e9 # typical steel value
+        Ip = klss*len_shaft/G
+        starr = np.zeros((2,19))
+        starr[:,0] = [0, len_shaft] # st file is always normalised
+        starr[:,1] = 0.5 # very light but not too light
+        starr[:,2] = 0.0 # no cg offset
+        starr[:,3] = 0.0 # no cg offset
+        # radius of gyration: very light but not too light
+        starr[:,4] = 0.001 #ri_x
+        starr[:,5] = 0.001 #ri_y
+        # shear center
+        starr[:,6] = 0.0 # no shear center offset
+        starr[:,7] = 0.0 # no shear center offset
+        # E, G: choose
+        starr[:,8] = 1e18
+        starr[:,9] = G
+        starr[:,10] = Ip/2 # area inertia's
+        starr[:,11] = Ip/2
+        starr[:,12] = Ip
+        # shear factor
+        starr[:,13] = 1000 # stiff in shear
+        starr[:,14] = 1000
+        starr[:,15] = 10 # cross-section area
+        starr[:,16] = 0.0 # structural pitch
+        starr[:,17] = 0.0 # elastic axis
+        starr[:,18] = 0.0
+        st = StFile()
+        st.main_data_sets = {1:{1:starr}}
+
+        # FIXME: INERTIATOCLUTCH??
+        inertiatoclutch = self.get_key('DTRAIN', 'INERTIATOCLUTCH')
+
+        # FIXME: what is the location of these inertia's?
+        # FIXME: what is the mass? It seems the mass is put elsewhere?
+
+        # cmass_gbx = {'Izz':Izz}
+        # cmass_gen = {'x':x, 'y':y, 'z':z, 'm':nacmas,
+        #              'Ixx':Ixx, 'Iyy':Iyy, 'Izz':Izz}
+
+        cmass_lss = {'x':0, 'y':0, 'z':0, 'm':0,'Ixx':0, 'Iyy':0, 'Izz':ginert}
+        # convert HSS inertia to LSS
+        cmass_hss = {'x':0, 'y':0, 'z':0, 'm':0,'Ixx':0, 'Iyy':0,
+                     'Izz':gbxinert*gratio*gratio}
+
+        return cmass_lss, cmass_hss, st
+
+    def get_nacelle(self):
+
+        # also get the nacell concentrated mass
+        nacmas = self.get_key('NMASS', 'NACMAS')[0,0]
+        # includes the nacelle structure and all the machinery within it.
+        # It does not include the rotor blades and hub. If the Direct Drive
+        # (see 4.1) option is selected, the mass of the generator (see 4.2)
+        # is also excluded.
+
+        # position of cg relative to tower top, nacelle mass coord system
+        # based on figure 4-7 of bladed manual 4.8
+        # FIXME: is this correct? z means SS, and at the edge of the nacelle?
+        nmx = self.get_key('NMASS', 'NMX')[0,0] # width?
+        nmy = self.get_key('NMASS', 'NMY')[0,0] # length?
+        nmz = self.get_key('NMASS', 'NMZ')[0,0] # height?
+        iyaw = self.get_key('NMASS', 'IYAW')[0,0]
+        inod = self.get_key('NMASS', 'INOD')[0,0]
+        iroll = self.get_key('NMASS', 'IROLL')[0,0]
+
+        # MSTART NGEOM
+        # NACW	 x.000
+        # NACL	 x.000
+        # NACH	 x.000
+        # NACZ	 0
+        # NACCD	 x.0
+        # NAERO	N
+
+        # height of the nacelle/towertop: the BLADED model does not have to
+        # have consistent shaft length with tilt angle arriving at the rotor
+        # center, so adjust towertop length in HAWC2 accordingly
+
+        # Tower height: from the ground or sea surface to the yaw bearing
+        # (only needed if the tower itself has not been defined)
+        towht = self.get_key('RCON', 'TOWHT')[0,0]
+
+        # Hub height: from the ground to the centre of the rotor
+        # (i.e. the intersection of the blade and shaft axes)
+        height = self.get_key('RCON', 'HEIGHT')[0,0]
+
+        # horizontal distance between rotor center and tower center line
+        ovrhng = self.get_key('RCON', 'OVRHNG')[0,0]
+        tilt = self.get_key('RCON', 'TILT')[0,0]
+
+        # FOR HAWC2, WE HAVE
+        len_shaft = ovrhng / np.cos(tilt)
+        len_towertop = height - towht - len_shaft*np.sin(tilt)
+        assert len_towertop > 0
+        # assert len_towertop + len_shaft*np.sin(tilt) + towht - height == 0
+
+        # concentrated mass in HAWC2
+        # 1  : nodenr
+        # 2-4: x, y, z offset
+        # 5  : mass
+        # 6-8: Ixx, Iyy, Izz
+        # 9-11: Ixy, Ixz, Iyz (optional)
+        # in HAWC2 we would be in the tower top coordinate system, and is
+        # usually the same as the global coordinate system
+        x = nmx
+        y = -nmy # move upwind wrt tower
+        z = -nmz # move upward
+        Ixx = inod
+        Iyy = iroll
+        Izz = iyaw
+
+        cmass_nacelle = {'x':x, 'y':y, 'z':z, 'm':nacmas,
+                         'Ixx':Ixx, 'Iyy':Iyy, 'Izz':Izz}
+
+        return cmass_nacelle, len_shaft, len_towertop
+
+    def get_towershadow_drag(self):
+        # TODO: tower shadow and tower aerodrag
+
+        # MSTART TOWSDW
+        # TSMODEL	x
+        # TDCORR	 x
+        # MAXDEF	 x.x
+        # WIDTH	 x
+        # LREF	 x.x
+
+        return
+
+    def get_control(self):
+
+        # MSTART CONTROL
+
+        # MSTART PCOEFF
+        # TSRMIN	 x
+        # TSRMAX	 x
+        # TSRSTP	 x.x
+        # PITCH	-x.0E-00
+        # PITCH_END	-x.0E-00
+        # PITCH_STEP	 0
+        # OMEGA	 x.654321
+
+        # MSTART IDLING
+
+        # AEROINFO
+
+        # pitch actuator rate limits
+        pitch_actuator = self.xmlroot.PitchActuator
+        if pitch_actuator.SetpointTrajectory.HasRateLimits:
+            # pitch_actuator.SetpointTrajectory.LowerRateLimit
+            # pitch_actuator.SetpointTrajectory.UpperRateLimit
+            print(pitch_actuator.SetpointTrajectory.UpperRateLimit*30/np.pi)
+
+        if pitch_actuator.SetpointTrajectory.HasAccelerationLimits:
+            # pitch_actuator.SetpointTrajectory.LowerAccelerationLimit
+            # pitch_actuator.SetpointTrajectory.UpperAccelerationLimit
+            print(pitch_actuator.SetpointTrajectory.UpperAccelerationLimit*30/np.pi)
+
+        # pitch_actuator.PositionResponse.values() # 1st or 2nd order
+        # pitch_actuator.PositionResponse.Frequency
+        # pitch_actuator.PositionResponse.Damping
+
+        # pitch_actuator.RateResponse.values() # 1st or 2nd order
+        # pitch_actuator.RateResponse.LagTime
+        # pitch_actuator.PositionResponse.Damping
+
+        # and more: bearing friction, etc
+
+        return
+
+    def convert(self, fname_htc):
+
+        fname_prj = os.path.basename(fname_htc).replace('.htc', '')
+
+        # for convenience, we start from an htc template that otherwise has the
+        # right structure, but has only 2 nodes for each body.
+        # assume it is in the same folder as where the HTC will be written to
+        basepath = os.path.dirname(fname_htc)
+        htc = HTCFile(pjoin(basepath, 'template.htc'))
+
+        # prj = self.PRJ_to_HAWC(pjoin(version, fname_prj))
+        # plot_bladed(prj)
+
+        # cone = prj.get_key('RCON', 'CONE')
+        # Lateral offset: the horizontal offset between the shaft and tower axes.
+
+        # Brake
+        # steps = prj.get_key('BRAKE', 'STEPS')
+        # TORQUE = prj.get_key('BRAKE', 'TORQUE')
+        # TIME = prj.get_key('BRAKE', 'TIME')
+
+        # -------------------------------------------------------------------------
+        # BLADES
+        # extract blade geometry, and add to c2_def section in the htc file
+        c2_blade = self.get_blade_c2_def()
+        htc = add_c2_arr_htc(htc, c2_blade, 'blade1')
+
+        # -------------------------------------------------------------------------
+        # NACELLE
+        cmass_nacelle, len_shaft, len_towertop = self.get_nacelle()
+        # set towertop length
+        # nacelle = htc['new_htc_structure'].get_subsection_by_name('towertop')
+        c2_tt = np.zeros((2,4))
+        c2_tt[1,2] = -len_towertop
+        htc = add_c2_arr_htc(htc, c2_tt, 'towertop')
+
+        # -------------------------------------------------------------------------
+        # attach nacelle mass/inertia to towertop, first node (yaw bearing)
+        towertop = htc['new_htc_structure'].get_subsection_by_name('towertop')
+        values = [v for k,v in cmass_nacelle.items()]
+        cmass = towertop.add_line('concentrated_mass', [1] + values)
+        key = cmass.location().split('/')[-1]
+        towertop[key].comments = 'nacelle mass (inc rotating) and inertia (non-rotating)'
+
+        # -------------------------------------------------------------------------
+        # set shaft length
+        c2_s = np.zeros((2,4))
+        c2_s[1,2] = len_shaft
+        htc = add_c2_arr_htc(htc, c2_s, 'shaft')
+
+        # -------------------------------------------------------------------------
+        # generator and gearbox inertia's, shaft torsional flexibility
+        cmass_lss, cmass_hss, st_shaft = self.get_drivertain(len_shaft)
+        st_shaft.save(pjoin(basepath, f'{fname_prj}_shaft.st'))
+
+        shaft = htc['new_htc_structure'].get_subsection_by_name('shaft')
+        # LSS inertia
+        values = [v for k,v in cmass_lss.items()]
+        cmass = shaft.add_line('concentrated_mass', [2] + values)
+        key = cmass.location().split('/')[-1]
+        shaft[key].comments = 'inertia LSS'
+        # Inertia of the HSS is already expressed wrt LSS
+        values = [v for k,v in cmass_hss.items()]
+        cmass = shaft.add_line('concentrated_mass', [2] + values)
+        gratio = self.get_key('DTRAIN', 'GRATIO')[0,0]
+        key = cmass.location().split('/')[-1]
+        shaft[key].comments = f'inertia HSS, expressed in LSS, GBR={gratio}'
+
+        # -------------------------------------------------------------------------
+        # set tilt angle
+        tilt = self.get_key('RCON', 'TILT')[0,0]
+        ori = htc['new_htc_structure'].orientation
+        # tilt is on the 2nd relative block, and set angle on 2nd set of eulerang
+        # check comments of the htc file:
+        # print(ori['relative__2']['mbdy2_eulerang__2'].comments)
+        ori['relative__2']['mbdy2_eulerang__2'].values = [tilt*180/np.pi, 0, 0]
+        ori['relative__2']['mbdy2_eulerang__2'].comments = 'tilt angle'
+
+        # coning angle
+        cone = self.get_key('RCON', 'CONE')[0,0]
+        # coning is on 3, 4 and 5 (all hub orientations)
+        for k in range(3,6):
+            # print(ori[f'relative__{k}']['mbdy2_eulerang__3'].comments)
+            ori[f'relative__{k}']['mbdy2_eulerang__3'].values = [-cone*180/np.pi, 0, 0]
+            ori[f'relative__{k}']['mbdy2_eulerang__3'].comments = 'cone angle'
+
+        # TODO: blade mounting sweep and cone as well??
+
+        # -------------------------------------------------------------------------
+        # HUB mass as concentrated mass on the shaft end
+        # FIXME: make sure we can get the hub root loads including these inertia's
+        cmass_hub, hublen = self.get_hub()
+        values = [v for k,v in cmass_hub.items()]
+        cmass = shaft.add_line('concentrated_mass', [2] + values)
+        key = cmass.location().split('/')[-1]
+        shaft[key].comments = 'Hub inertia and mass'
+
+        # -------------------------------------------------------------------------
+        # hub lenght
+        c2_hub = np.zeros((2,4))
+        c2_hub[1,2] = hublen
+        htc = add_c2_arr_htc(htc, c2_hub, 'hub1')
+
+        # -------------------------------------------------------------------------
+        # TOWER ST FILE and C2_DEF section
+        # WARNING: only includes nodes that are at centerline, side-elements are
+        # ignored
+        # TODO: compare between TGEOM/TMASS and TSTIFF approaches, it seems that
+        # the shear stiffness is not the same, see check_tower()
+        c2_tow, st = self.get_tower_c2_def_st()
+        st.save(pjoin(basepath, f'{fname_prj}_tower.st'))
+        htc = add_c2_arr_htc(htc, c2_tow, 'tower')
+
+        # hub height
+        htc.wind.center_pos0.values = [0, 0, -self.get_key('RCON', 'HEIGHT')[0,0]]
+
+        # -------------------------------------------------------------------------
+        # set all file names correct
+        tower = htc['new_htc_structure'].get_subsection_by_name('tower')
+        tower.timoschenko_input.filename.values[0] = f'data/{fname_prj}_tower.st'
+
+        shaft = htc['new_htc_structure'].get_subsection_by_name('shaft')
+        shaft.timoschenko_input.filename.values[0] = f'data/{fname_prj}_shaft.st'
+
+        blade1 = htc['new_htc_structure'].get_subsection_by_name('blade1')
+        blade1.timoschenko_input.filename.values[0] = f'data/{fname_prj}_blade.st'
+
+        for body_name in ['towertop', 'hub1']:
+            body = htc['new_htc_structure'].get_subsection_by_name(body_name)
+            body.timoschenko_input.filename.values[0] = 'data/template.st'
+
+        htc['aero'].ae_filename.values[0] = f'data/{fname_prj}.ae'
+        htc['aero'].pc_filename.values[0] = f'data/{fname_prj}.pc'
+
+        # -------------------------------------------------------------------------
+        # SAVE HTC FILE
+        htc.save(pjoin(basepath, f'{fname_htc}'))
+
+        # -------------------------------------------------------------------------
+        # other data files
+
+        # Blade st file
+        st = self.get_blade_st()
+        st.save(pjoin(basepath, f'{fname_prj}_blade.st'))
+
+        # extract profile coefficients and save to pc file
+        pc = self.get_pc(setname_filter='cleanVG')
+        pc.save(pjoin(basepath, f'{fname_prj}.pc'))
+
+        # extract aerodynamic layout, and safe to ae file
+        ae = self.get_ae()
+        ae.save(pjoin(basepath, f'{fname_prj}.ae'))
+
+        # # BLADE c2_def
+        # c2_def = htc['new_htc_structure'].get_subsection_by_name('blade1')['c2_def']
+        # # iterate over all nsec lines, but ignore the line that gives the number of points
+        # for i, sec in enumerate(c2_def.keys()[1:]):
+        #     c2_arr[i,:] = c2_def[sec].values[1:]
+
+        # curvedlength = np.zeros(len(c2_arr[:,2]))
+        # for i in range(1,len(c2_arr[:,2])):
+        #     curvedlength[i]=curvedlength[i-1]+np.linalg.norm(c2_arr[i,0:3]-c2_arr[i-1,0:3])
+
+        # c2_arr = np.zeros((len(c2_def.keys())-1,4))
+        # # iterate over all nsec lines, but ignore the line that gives the number of points
+        # for i, sec in enumerate(c2_def.keys()[1:]):
+        #     c2_arr[i,:] = c2_def[sec].values[1:]
--- a/wetb/bladed/readprj.py
+++ b/wetb/bladed/readprj.py
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Wed Jun 23 11:07:50 2021
+
+@author: dave
+"""
+
+# from os.path import join as pjoin
+
+import numpy as np
+from lxml import objectify#, etree)
+
+# from matplotlib import pyplot as plt
+# import pandas as pd
+
+# import wetb
+from wetb.prepost import misc
+# from wetb.hawc2 import (HTCFile, AEFile, PCFile, StFile)
+
+
+# TODO: this dictionary utility could go to a dict sub-class in misc?
+def find_next_unique_key(d, key):
+    k = 1
+    while True:
+        k += 1
+        if key + f'__{k}' not in d:
+            key = key + f'__{k}'
+            return key
+
+
+# TODO: this dictionary utility could go to a dict sub-class in misc?
+def add_new_unique_key(key, values, d):
+    """
+    Add key:values to dictionary d. If key already occurs, append __x
+    where x is the number of previous occurances of key(__x).
+
+    Parameters
+    ----------
+    key : str
+        DESCRIPTION.
+    values : list
+        DESCRIPTION.
+    d : dict
+        Dictionary.
+
+    Returns
+    -------
+    key : str
+        DESCRIPTION.
+    d : dict
+        DESCRIPTION.
+
+    """
+    if key in d:
+        key = find_next_unique_key(d, key)
+    d[key] = [values]
+    return key, d
+
+
+class ReadBladedProject:
+
+    def __init__(self, fname):
+
+        with open(fname, encoding='utf-8') as fobj:
+            xml_str = fobj.read().encode('utf-8')
+        self.bd, self.xmlroot = self.parse_bladeddata(xml_str)
+
+        self.set_attr_and_check()
+
+        # some things are just a little different
+        # TMASS has a list of materials and their properties
+        # get rid of the quotes
+        # tmp = [el.replace("'", '') for el in self.bd['TMASS']['MATERIAL'][0]]
+        # self.bd['TMASS']['MATERIAL'] = tmp
+        unique_mat = set(self.get_key('TMASS', 'MATERIAL').flatten().tolist())
+        self.tow_mat_prop = {k:self.get_key('TMASS', k) for k in unique_mat}
+        # material_props = {}
+        # for k in unique_mat:
+        #     material_props[k.replace("'", '')] = self.get_key('TMASS', k)
+
+    def parse_bladeddata(self, xml_str):
+        """
+        The XML field BladedData contains what seems like the main core input
+        data for BLADED, and formatted in some structured way.
+
+        Parameters
+        ----------
+        xml_str : TYPE
+            DESCRIPTION.
+
+        Returns
+        -------
+        bd : TYPE
+            DESCRIPTION.
+        xmlroot : TYPE
+            DESCRIPTION.
+
+        """
+
+        # root = etree.fromstring(xml)
+        # elems = root.getchildren()
+        # bladeddata = elems[1].text
+
+        xmlroot = objectify.fromstring(xml_str)
+
+        # the non-xml formatted BLADED model is situated in a non-xml field
+        # called BladedData
+        bd = {}
+        mstart = None
+        for i, line in enumerate(xmlroot.BladedData.text.split('\n')):
+
+            # TODO: values embedded in double quotes (") can contain entire
+            # MSTART/MEND sub-sections (so embedded sections)
+
+            # split replace tabs with spaces, split on spaces, remove empty
+            linels = misc.remove_items(line.replace('\t', ' ').split(' '), '')
+            # commas can also be separators, in addition to spaces
+            linels2 = []
+            for k in linels:
+                linels2.extend(k.split(','))
+            linels = misc.remove_items(linels2, '')
+
+            # ignore empty lines
+            if len(linels) < 1:
+                continue
+
+            # entries can be numbers if the previous key contains multiple data points
+            try:
+                float(linels[0])
+                el0isnum = True
+            except ValueError:
+                el0isnum = False
+
+            # start of a sub-section that contains (non-unique) keys as well
+            if linels[0].upper().startswith('MSTART'):
+                mtag = linels[-1]
+                mstart = {}
+
+            # at the end of the sub-section add the sub-section to the main dict
+            elif linels[0].upper().startswith('MEND'):
+                # FIXME: implement MSTART sections embedded in double quoted values
+                try:
+                    # if the section key is not unique, make it so by appending __X
+                    if mtag in bd:
+                        mtag = find_next_unique_key(bd, mtag)
+                    bd[mtag] = mstart
+                except UnboundLocalError:
+                    print('warning: ignored embedded mstart/mend section')
+                    print(f'at line: {i+1}')
+                mstart = None
+
+            # if we are under a sub-section
+            elif mstart is not None:
+                # if the line contains a keyword
+                if not el0isnum:
+                    tag, mstart = add_new_unique_key(linels[0], linels[1:], mstart)
+                # line is datapoint that needs to be added to key that occured before
+                else:
+                    mstart[tag].append(linels)
+
+            # add numerical values to key that occured before
+            elif el0isnum:
+                 bd[tag].append(linels)
+
+            else:
+                tag, bd = add_new_unique_key(linels[0], linels[1:], bd)
+
+        return bd, xmlroot
+
+    def get_key(self, key1, key2=False):
+        """
+        Get key from the BladedData CDATA section and format to int32 or
+        float32 numpy arrays if possible withouth precision loss.
+
+        Parameters
+        ----------
+        key1 : str
+            DESCRIPTION.
+        key2 : str, optional
+            DESCRIPTION. The default is False.
+
+
+        Returns
+        -------
+        numpy.array
+            Values from key1/key2 formatted as a numpy array. Converted to
+            numpy.int32, numpy.float32 if possible withouth precision loss,
+            otherwise an object array is returned.
+
+        """
+
+        if key1 not in self.bd:
+            raise KeyError(f'{key1} not found in BLADED file.')
+
+        if key2 is not False:
+            if key2 not in self.bd[key1]:
+                raise KeyError(f'{key2} not found in MSTART {key1} of BLADED file.')
+            data = self.bd[key1][key2]
+        else:
+            data = self.bd[key1]
+
+        # in case we defined a mstart key
+        if isinstance(data, dict):
+            return data
+
+        # i ,j = len(data), len(data[0])
+
+        # FIXME: this is a very expensive way of converting it, but it might
+        # not matter at all since very little model data is actually considered
+        data_arr = np.array(data)
+        try:
+            data_arr = data_arr.astype(np.int32)
+        except ValueError:
+            try:
+                data_arr = data_arr.astype(np.float32)
+            except ValueError:
+                pass
+
+        return data_arr
+
+        # return np.array(data, dtype=np.float32)
+
+        # if isinstance(data[0], list) and len(data[0]) == 1:
+        #     data = float(data[0])
+        # if isinstance(data[0], list) and len(data[0]) > 1:
+        #     data_arr = np.array(data, dtype=np.float32)
+
+    def set_attr_and_check(self):
+        """Check that BGEOMMB, BMASSMB, BSTIFFMB has indeed the same node
+        repated twice every time.
+        """
+
+        # self.bd['BGEOMMB'].keys(), but only those relevant
+        keysg = ['RJ', 'DIST', 'REF_X', 'REF_Y', 'CHORD', 'TWIST', 'CE_X',
+                 'CE_Y', 'BTHICK', 'FOIL', 'MOVING']
+        nbe = self.get_key('BGEOMMB', 'NBE')[0,0]
+
+        # self.bd['BMASSMB'].keys(), but only those relevant
+        keysm = ['CM_X', 'CM_Y', 'MASS', 'SINER', 'RGRATIO', 'BETA_M']
+
+        # self.bd['BSTIFFMB'].keys(), but only those relevant
+        keyss = ['EIFLAP', 'EIEDGE', 'BETA_S', 'GJ', 'CS_X', 'CS_Y', 'GAFLAP',
+                 'GAEDGE']
+
+        mkeys = ['BGEOMMB', 'BMASSMB', 'BSTIFFMB']
+        for mkey, keys in zip(mkeys, [keysg, keysm, keyss]):
+            for key in keys:
+                res = self.get_key(mkey, key)
+                try:
+                    assert np.allclose(res[0,0::2], res[0,1::2])
+                except TypeError:
+                    # allclose doesn't make sense for text arrays
+                    assert np.compare_chararrays(res[0,0::2], res[0,1::2],
+                                                 '==', True).all()
+                assert res.shape[1]==nbe
+                if hasattr(self, key.lower()):
+                    raise UserWarning(key, 'already exists')
+                setattr(self, key.lower(), res[0,0::2])
+
+    def print_xml_tree(self, fname):
+        """For discovery purposes: print full tree + values/text
+
+
+        Parameters
+        ----------
+        fname : TYPE
+            DESCRIPTION.
+
+        Returns
+        -------
+        None.
+
+        """
+
+        # def print_rec(root):
+        #     # if hasattr(root, 'getparent'):
+        #     #     print(root.getparent().tag.title, end='.')
+        #     #     print_rec(root.getparent())
+        #     for el in root.getchildren():
+        #         print(print_rec(el))
+
+        # def print_root_tree(root):
+        #     root.getparent()
+        #     print()
+
+        # tree = etree.fromstring(xml_str)
+        # els = tree.xpath('/')
+        # for el in els:
+        #     print(el)
+
+        # tree = etree.fromstring(xml_str)
+
+        # xmlroot = objectify.fromstring(xml_str)
+        # with open(fname+'.structure', 'w') as f:
+        #     for line in xmlroot.descendantpaths():
+        #         f.write(line + '\n')
+
+        # Recursive XML parsing python using ElementTree
+        # https://stackoverflow.com/q/28194703/3156685
+
+        roottree = self.xmlroot.getroottree()
+        def print_tree_recursive(root):
+            # print(roottree.getpath(root), end=' : ')
+            # print(root.tag.title())
+            f.write(roottree.getpath(root) + ' : ')
+            f.writelines(root.values())
+            if root.text is not None:
+                f.write(' : ' + root.text)
+            f.write('\n')
+            for elem in root.getchildren():
+                print_tree_recursive(elem)
+        with open(fname+'.structure.value', 'w') as f:
+            for el in self.xmlroot.getchildren():
+                if el.tag.title()!='Bladeddata':
+                    print_tree_recursive(el)
+
+    # def get_frequencies(self):
+    #     """
+    #     """
+
+    #     blades = self.xmlroot.NewGUI.Turbine.Blades.FlexibilityModel
+
+    #     blades.Settings.ModesWithDampingDefined
+    #     blades.PartModes # with lots of BladeModeContainer
+    #     blades.WholeBladeModes
+    #     blades.WholeBladeModes.WholeBladeMode
+
+    #     blades.WholeBladeModes.WholeBladeMode.Name
+    #     blades.WholeBladeModes.WholeBladeMode.Frequency
+    #     blades.WholeBladeModes.WholeBladeMode.Damping
+    #     blades.WholeBladeModes.WholeBladeMode.Components
+    #     len(blades.WholeBladeModes.WholeBladeMode.Components.getchildren())
+    #     # 60 elements
+
+    #     for wholeblademode in blades.WholeBladeModes.iterchildren():
+    #         print(wholeblademode.Name)
+    #         print(wholeblademode.Frequency, wholeblademode.Damping)
+
+    #     tower = self.xmlroot.NewGUI.Turbine.SupportStructure.FlexibilityModel
+    #     for towermode in tower.Modes.getchildren():
+    #         print(towermode.Description)
+    #         towermode.Frequency
+    #         towermode.Damping
--- a/wetb/bladed/template.htc
+++ b/wetb/bladed/template.htc
+;
+begin simulation;
+  time_stop    1;
+  solvertype   1 ;    (newmark)
+  on_no_convergence continue ;
+  convergence_limits 1E3 1.0 1E-7 ;
+  logfile case0.log ;
+  begin newmark;
+    deltat    0.02;
+  end newmark;
+end simulation;
+;----------------------------------------------------------------------------------------------------------------------
+begin new_htc_structure;
+;----------------------------------------------------------------------------------------------------------------------
+  begin main_body;
+    name        tower ;
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution     c2_def ;
+    damping_posdef   0.0 0.0 0.0 0.0E-00 0.0E-00 0.0E-00;
+     begin timoschenko_input;
+      filename template.st;
+      set 1 1 ;
+    end timoschenko_input;
+    begin c2_def;
+      nsec 2;
+      sec   1   0   0     0.00  0 ;
+      sec   2   0   0  -100.00  0 ;
+     end c2_def ;
+    end main_body;
+;
+  begin main_body;
+    name        towertop ;
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution     c2_def ;
+    damping_posdef   0.0 0.0 0.0 0.0E-00 0.0E-00 0.0E-00;
+    begin timoschenko_input;
+      filename template.st;
+      set 1 1 ;
+    end timoschenko_input;
+    begin c2_def;
+      nsec 2;
+      sec 1 0.0 0.0  0.0    0.0 ;
+      sec 2 0.0 0.0 -1.00   0.0 ;
+    end c2_def ;
+  end main_body;
+;
+  begin main_body;
+    name        shaft ;
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution     c2_def ;
+    damping_posdef   0.0 0.0 0.0 0.0E-00 0.0E-00 0.0E-00;
+    begin timoschenko_input;
+      filename template.st;
+      set 1 1 ;
+    end timoschenko_input;
+    begin c2_def;
+      nsec 2;
+      sec 1 0.0 0.0 0.0 0.0 ;
+      sec 2 0.0 0.0 1.0 0.0 ;
+    end c2_def ;
+  end main_body;
+;
+  begin main_body;
+    name        hub1 ;
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution     c2_def ;
+    damping_posdef   0.0 0.0 0.0 0.0E-00 0.0E-00 0.0E-00;
+    begin timoschenko_input;
+      filename template.st;
+      set 1 1 ;
+    end timoschenko_input;
+    begin c2_def;
+      nsec 2;
+      sec 1 0.0 0.0 0.0 0.0 ;
+      sec 2 0.0 0.0 1.0 0.0 ;
+    end c2_def ;
+  end main_body;
+;
+  begin main_body;
+    name           hub2 ;
+    copy_main_body hub1;
+  end main_body;
+;
+  begin main_body;
+    name           hub3 ;
+    copy_main_body hub1 ;
+  end main_body;
+;
+  begin main_body;
+    name        blade1 ;
+    type        timoschenko ;
+    nbodies     1 ;
+    node_distribution    c2_def;
+    damping_posdef   0.0 0.0 0.0 0.0E-00 0.0E-00 0.0E-00;
+    begin timoschenko_input ;
+      filename template.st;
+      set 1 1 ;
+    end timoschenko_input;
+    begin c2_def;
+        nsec 2;
+        sec 1    0.00000E+00     0.00000E+00    0.00000E+00   0.00000E+00 ;
+        sec 2    0.00000E+00     0.00000E+00    8.00000E+01   0.00000E+00 ;
+     end c2_def ;
+   end main_body;
+;
+  begin main_body;
+    name           blade2 ;
+    copy_main_body blade1;
+  end main_body;
+;
+  begin main_body;
+    name           blade3 ;
+    copy_main_body blade1 ;
+  end main_body;
+;----------------------------------------------------------------------------------------------------------------------
+;
+  begin orientation;
+    begin base;
+      mbdy   tower;
+      inipos        0.0 0.0 0.0 ;
+      body_eulerang 0.0 0.0 0.0;
+    end base;
+    begin relative;
+      mbdy1  tower last;
+      mbdy2  towertop 1;
+      mbdy2_eulerang 0.0 0.0 0.0;
+    end relative;
+    begin relative;
+      mbdy1  towertop last;
+      mbdy2  shaft 1;
+      mbdy2_eulerang 90.0 0.0 0.0;
+      mbdy2_eulerang 0.0 0.0 0.0; [tilt] 0.0 0.0
+      mbdy2_eulerang 0.0 0.0 0.0; 0.0 0.0 [azi]
+      mbdy2_ini_rotvec_d1 0.0 0.0 -1.0 0.0 ; [init_wr];
+    end relative;
+    begin relative;
+      mbdy1  shaft last;
+      mbdy2  hub1 1;
+      mbdy2_eulerang -90.0 0.0 0.0;
+      mbdy2_eulerang 0.0 180.0 0.0;
+      mbdy2_eulerang 2.5 0.0 0.0; [cone] 0.0 0.0
+    end relative;
+    begin relative;
+      mbdy1  shaft last;
+      mbdy2  hub2 1;
+      mbdy2_eulerang -90.0 0.0 0.0;
+      mbdy2_eulerang 0.0 60.0 0.0;
+      mbdy2_eulerang 2.5 0.0 0.0;
+    end relative;
+    begin relative;
+      mbdy1  shaft last;
+      mbdy2  hub3 1;
+      mbdy2_eulerang -90.0 0.0 0.0;
+      mbdy2_eulerang 0.0 -60.0 0.0;
+      mbdy2_eulerang 2.5 0.0 0.0;
+    end relative;
+    begin relative;
+      mbdy1  hub1 last;
+      mbdy2  blade1 1;
+      mbdy2_eulerang 0.0 0.0 0;
+    end relative;
+    begin relative;
+      mbdy1  hub2 last;
+      mbdy2  blade2 1;
+      mbdy2_eulerang 0.0 0.0 0.0;
+    end relative;
+    begin relative;
+      mbdy1  hub3 last;
+      mbdy2  blade3 1;
+      mbdy2_eulerang 0.0 0.0 0.0;
+    end relative;
+   end orientation;
+;----------------------------------------------------------------------------------------------------------------------
+  begin constraint;
+    begin fix0;
+      body tower;
+    end fix0;
+    begin fix1;
+       mbdy1 tower last ;
+       mbdy2 towertop 1;
+    end fix1;
+    begin bearing1;
+      name  shaft_rot;
+      mbdy1 towertop last;
+      mbdy2 shaft 1;
+      bearing_vector 2 0.0 0.0 -1.0; (0=global.1=mbdy1.2=mbdy2)
+    end bearing1;
+;     begin bearing3;
+;       name  shaft_rot;
+;       mbdy1 towertop last;
+;       mbdy2 shaft 1;
+;       bearing_vector 2 0.0 0.0 -1.0; (0=global.1=mbdy1.2=mbdy2)
+;       omegas 0.0 ;
+;    end bearing3;
+    begin fix1;
+       mbdy1 shaft last ;
+       mbdy2 hub1 1;
+    end fix1;
+    begin fix1;
+       mbdy1 shaft last ;
+       mbdy2 hub2 1;
+    end fix1;
+    begin fix1;
+       mbdy1 shaft last ;
+       mbdy2 hub3 1;
+    end fix1;
+    begin bearing2;
+       name pitch1;
+       mbdy1 hub1 last;
+       mbdy2 blade1 1;
+       bearing_vector 2 0.0 0.0 -1.0;
+    end bearing2;
+    begin bearing2;
+       name pitch2;
+       mbdy1 hub2 last;
+       mbdy2 blade2 1;
+       bearing_vector 2 0.0 0.0 -1.0;
+    end bearing2;
+    begin bearing2;
+       name pitch3;
+       mbdy1 hub3 last;
+       mbdy2 blade3 1;
+       bearing_vector 2 0.0 0.0 -1.0;
+    end bearing2;
+  end constraint;
+end new_htc_structure;
+;----------------------------------------------------------------------------------------------------------------------
+begin wind ;
+  density                 1.225 ;
+  wsp                     1;
+  tint                    0.01;
+  horizontal_input        1     ;
+  windfield_rotations     0.0 0.0 0.0 ;    yaw, tilt, rotation
+  center_pos0             0.0 0.0 -119 ; hub heigth
+  shear_format            3 0.2;
+  turb_format             0     ;  0=none, 1=mann,2=flex
+  tower_shadow_method     3     ;  0=none, 1=potential flow, 2=jet
+;  scale_time_start       100;
+;  wind_ramp_factor   0.0 50 0.0 1.0 ;
+;
+  begin tower_shadow_potential_2;
+    tower_mbdy_link tower;
+    nsec  2;
+    radius      0.0  4.15 ;
+    radius     115.63 2.75 ;
+  end tower_shadow_potential_2;
+end wind;
+;
+begin aerodrag ;
+  begin aerodrag_element ;
+    mbdy_name tower;
+    aerodrag_sections uniform 10 ;
+    nsec 2 ;
+    sec   0.0 0.6 1.0 ;  tower bottom
+    sec 100.0 0.6 1.0 ;  tower top
+  end aerodrag_element;
+;
+  begin aerodrag_element ; Nacelle drag side
+    mbdy_name shaft;
+    aerodrag_sections uniform 2 ;
+    nsec 2 ;
+    sec 0.0   0.8 10.0 ;
+    sec 1.00  0.8 10.0 ;
+  end aerodrag_element;
+end aerodrag;
+;
+begin aero ;
+  nblades  3;
+  hub_vec shaft -3 ;
+  link 1 mbdy_c2_def blade1;
+  link 2 mbdy_c2_def blade2;
+  link 3 mbdy_c2_def blade3;
+  ae_filename        template.ae ;
+  pc_filename        template.pc ;
+  induction_method   0 ;     0=none, 1=normal
+  aerocalc_method    1 ;     0=ingen aerodynamic, 1=med aerodynamic
+  aerosections       2;
+  ae_sets            1 1 1;
+  tiploss_method     1 ;     0=none, 1=prandtl
+  dynstall_method    2 ;     0=none, 1=stig øye method,2=mhh method
+end aero ;
+;----------------------------------------------------------------------------------------------------------------------
+begin output;
+  filename case0;
+  time 0 100;
+  data_format gtsdf;
+  buffer 9999;
+;
+  general time;
+  constraint bearing1 shaft_rot 2; angle and angle velocity
+  constraint bearing2 pitch1 5;    angle and angle velocity
+  constraint bearing2 pitch2 5;    angle and angle velocity
+  constraint bearing2 pitch3 5;    angle and angle velocity
+  aero omega ;
+  aero torque;
+  aero power;
+  aero thrust;
+  wind free_wind 1 0.0 0.0 -100;
+  ; Moments:
+  mbdy momentvec tower  1 1  tower # tower base ;
+  mbdy momentvec tower  1 2 tower # tower yaw bearing ;
+  mbdy momentvec shaft  1 1  shaft # main bearing ;
+  mbdy momentvec blade1 1 1 blade1 # blade 1 root ;
+  mbdy momentvec blade2 1 1 blade2 # blade 2 root ;
+  mbdy momentvec blade3 1 1 blade3 # blade 3 root ;
+  mbdy momentvec blade1 1 2 local # blade 1 50% local e coo ;
+  mbdy momentvec blade2 1 2 local # blade 2 50% local e coo ;
+  mbdy momentvec blade3 1 2 local # blade 3 50% local e coo ;
+  ; Displacements and accellerations
+  mbdy state pos tower 1 1.0 global only 1 # Tower top FA displ;
+  mbdy state pos tower 1 1.0 global only 2 # Tower top SS displ;
+  mbdy state acc tower 1 1.0 global only 1 # Tower top FA acc;
+  mbdy state acc tower 1 1.0 global only 2 # Tower top SS acc;
+  mbdy state pos blade1  1 1.0 global # gl blade 1 tip pos ;
+  mbdy state pos blade2  1 1.0 global # gl blade 2 tip pos ;
+  mbdy state pos blade3  1 1.0 global # gl blade 3 tip pos ;
+  mbdy state pos blade1  1 1.0 blade1 # blade 1 tip pos ;
+end output;
+;
+exit;
--- a/wetb/bladed/template.st
+++ b/wetb/bladed/template.st
+1  number of sets, Nset
+-----------------
+#1 Nset number 1
+================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================
+           r [0]                    m [1]                   x_cg [2]                y_cg [3]                    ri_x [4]                ri_y [5]                x_sh [6]                y_sh [7]                    E [8]                   G [9]                   I_x [10]                 I_y [11]                  K [12]                  k_x [13]                 k_y [14]                  A [15]                 pitch [16]                x_e [17]                 y_e [18]
+================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================
+$1 2 subset number 1
+   0.000000000000000e+00    0.500000000000000e+00    0.000000000000000e+00    0.000000000000000e+00    0.500000000000000e+00    0.500000000000000e+00    0.000000000000000e+00    0.000000000000000e+00    1.000000000000000e+14    1.000000000000000e+13    1.000000000000000e+00    1.000000000000000e+00    1.000000000000000e+00    1.000000000000000e+00    1.000000000000000e+00    1.000000000000000e+00    0.000000000000000e+00    0.000000000000000e+00    0.000000000000000e+00
+   1.000000000000000e+00    0.500000000000000e+00    0.000000000000000e+00    0.000000000000000e+00    0.500000000000000e+00    0.500000000000000e+00    0.000000000000000e+00    0.000000000000000e+00    1.000000000000000e+14    1.000000000000000e+13    1.000000000000000e+00    1.000000000000000e+00    1.000000000000000e+00    1.000000000000000e+00    1.000000000000000e+00    1.000000000000000e+00    0.000000000000000e+00    0.000000000000000e+00    0.000000000000000e+00
+
--- a/wetb/control/control.py
+++ b/wetb/control/control.py
@@ -5,10 +5,6 @@ Created on Thu Aug 04 09:24:51 2016
 @author: tlbl
 """

-from __future__ import division
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import absolute_import

 import numpy as np


--- a/wetb/control/tests/test_control.py
+++ b/wetb/control/tests/test_control.py
@@ -5,12 +5,6 @@ Created on Thu Aug 04 11:09:43 2016
 @author: tlbl
 """

-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import division
-from __future__ import absolute_import
-from future import standard_library
-standard_library.install_aliases()
 import unittest

 from wetb.control import control
No results found